{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 第四题：神经网络：三层感知机"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "实现内容：\n",
    "1. 实现一个三层感知机\n",
    "2. 对手写数字数据集进行分类\n",
    "3. 绘制损失值变化曲线\n",
    "4. 完成kaggle MNIST手写数字分类任务，根据给定的超参数训练模型，完成表格的填写"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "在这道题中，我们要实现一个三层感知机\n",
    "\n",
    "<img src=\"https://davidham3.github.io/blog/2018/09/11/logistic-regression/Fig2.png\" ,width=600>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 前向传播\n",
    "\n",
    "我们实现一个最简单的三层感知机，一个输入层，一个隐藏层，一个输出层，隐藏层单元个数为$h$个，输出层有$K$个单元。\n",
    "\n",
    "1. 我们将第一层的输入，定义为$X \\in \\mathbb{R}^{n \\times m}$，n个样本，m个特征。  \n",
    "2. 输入层到隐藏层之间的权重(weight)与偏置(bias)，分别为$W_1 \\in \\mathbb{R}^{m \\times h}$，$b_1 \\in \\mathbb{R}^{1 \\times h}$。  \n",
    "3. 隐藏层到输出层的权重和偏置分为别$W_2 \\in \\mathbb{R}^{h \\times K}$，$b_2 \\in \\mathbb{R}^{1 \\times K}$。\n",
    "\n",
    "隐藏层的激活函数选用ReLU\n",
    "\n",
    "$$\n",
    "\\mathrm{ReLU}(x) = \\max (0, x)\n",
    "$$\n",
    "\n",
    "我们用$H_1$表示第一个隐藏层的输出值，$O$表示输出层的输出值，这样，前向传播即可定义为\n",
    "\n",
    "$$\n",
    "Z = XW_1 + b_1\\\\\n",
    "H_1 = \\mathrm{ReLU}(Z)\\\\\n",
    "O = H_1 W_2 + b_2\n",
    "$$\n",
    "\n",
    "其中，$H_1 \\in \\mathbb{R}^{n \\times h}$，$O \\in \\mathbb{R}^{n \\times K}$。\n",
    "\n",
    "**注意：这里我们其实是做了广播，将$b_1$复制了$n-1$份后拼接成了维数为$n \\times h$的矩阵，同理，$b_2$也做了广播，拼成了$n \\times K$的矩阵。**\n",
    "\n",
    "最后一层的输出，使用softmax函数激活，得到神经网络计算出的各类的概率值：\n",
    "\n",
    "$$\n",
    "\\begin{aligned}\n",
    "\\hat{y_i} & = \\mathrm{softmax}(O_i)\\\\\n",
    "& = \\frac{\\exp{(O_i)}}{\\sum^{K}_{k=1} \\exp{(O_k)}}\n",
    "\\end{aligned}\n",
    "$$\n",
    "\n",
    "其中，$\\hat{y_i}$表示第$i$类的概率值，也就是输出层第$i$个神经元经$\\mathrm{softmax}$激活后的值。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 损失函数\n",
    "\n",
    "损失函数使用交叉熵损失函数：\n",
    "$$\\mathrm{cross\\_entropy}(y, \\hat{y}) = -\\sum^{K}_{k=1}y_k \\log{(\\hat{y_k})}$$\n",
    "\n",
    "这样，$n$个样本的平均损失为：\n",
    "$$\n",
    "\\mathrm{loss} = - \\frac{1}{n} \\sum_n \\sum^{K}_{k=1} y_k \\log{(\\hat{y_k})}\n",
    "$$\n",
    "\n",
    "**注意，这里我们的提到的$\\log$均为$\\ln$，在numpy中为**`np.log`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 反向传播\n",
    "\n",
    "我们使用梯度下降训练模型，求解方式就是求出损失函数对参数的偏导数，即参数的梯度，然后将参数减去梯度乘以学习率，进行参数的更新。\n",
    "$$\n",
    "W := W - \\alpha \\frac{\\partial \\mathrm{loss}}{\\partial W}\n",
    "$$\n",
    "其中，$\\alpha$是学习率。\n",
    "\n",
    "在这道题中，交叉熵损失函数的求导比较麻烦，我们先求神经网络的输出层的偏导数，写成链式法则的形式：\n",
    "\n",
    "$$\n",
    "\\frac{\\partial \\mathrm{loss}}{\\partial O_i} = \\frac{\\partial \\mathrm{loss}}{\\partial \\hat{y}} \\frac{\\partial \\hat{y}}{\\partial O_i}\n",
    "$$\n",
    "\n",
    "首先求解第一项：\n",
    "$$\n",
    "\\frac{\\partial \\mathrm{loss}}{\\partial \\hat{y}} = - \\frac{1}{n} \\sum_n \\sum^{K}_{k=1} y_k \\frac{1}{\\hat{y_k}}\n",
    "$$\n",
    "\n",
    "然后求解第二项，因为$\\hat{y_k}$的分母是$\\sum_k \\exp{(O_k)}$，里面包含$O_i$，所以每一个$\\hat{y_k}$的分母都包含$O_i$，这就要求反向传播的时候需要考虑这$K$项，将这$K$项的偏导数加在一起。\n",
    "\n",
    "这$K$项分别为：$\\frac{\\exp{(O_1)}}{\\sum_k \\exp{(O_k)}}$，$\\frac{\\exp{(O_2)}}{\\sum_k \\exp{(O_k)}}$，...，$\\frac{\\exp{(O_i)}}{\\sum_k \\exp{(O_k)}}$，...，$\\frac{\\exp{(O_k)}}{\\sum_k \\exp{(O_k)}}$。\n",
    "\n",
    "显然，这里只有分子带有$O_i$的这项与其他的项不同，因为分子和分母同时包含了$O_i$，而其他的项只有分母包含了$O_i$。\n",
    "\n",
    "这就需要在求解$\\frac{\\partial \\hat{y}}{\\partial O_i}$的时候分两种情况讨论\n",
    "1. 分子带$O_i$\n",
    "2. 分子不带$O_i$\n",
    "\n",
    "第一种情况，当分子含有$O_i$时：\n",
    "\n",
    "$$\n",
    "\\begin{aligned}\n",
    "\\frac{\\partial \\hat{y_i}}{\\partial O_i} & = \\frac{\\partial \\hat{y_i}}{\\partial O_i}\\\\\n",
    "& = \\frac{\\exp{(O_i)} (\\sum^{K}_{k=1} \\exp{(O_k)}) - (\\exp{(O_i)})^2 }{(\\sum^{K}_{k=1} \\exp{(O_k)})^2}\\\\\n",
    "& = \\frac{\\exp{(O_i)}}{\\sum^{K}_{k=1} \\exp{(O_k)}} \\frac{\\sum^{K}_{k=1} \\exp{(O_k)} - \\exp{(O_i)}}{\\sum^{K}_{k=1} \\exp{(O_k)}}\\\\\n",
    "& = \\hat{y_i} ( 1 - \\hat{y_i} )\n",
    "\\end{aligned}\n",
    "$$\n",
    "\n",
    "第二种情况，当分子不含$O_i$时，我们用$j$表示当前项的下标：\n",
    "\n",
    "$$\n",
    "\\begin{aligned}\n",
    "\\frac{\\partial \\hat{y_j}}{\\partial O_i} & = \\frac{- \\exp{(O_j)} \\exp{(O_i)}}{(\\sum^{K}_{k=1} \\exp{(O_k)})^2}\\\\\n",
    "& = - \\hat{y_j} \\hat{y_i}\n",
    "\\end{aligned}\n",
    "$$\n",
    "\n",
    "这样，$\\mathrm{loss}$对$O_i$的偏导数即为：\n",
    "$$\\begin{aligned}\n",
    "\\frac{\\partial \\mathrm{loss}}{\\partial O_i} & = \\frac{\\partial \\mathrm{loss}}{\\partial \\hat{y}} \\frac{\\partial \\hat{y}}{\\partial O_i}\\\\\n",
    "& = (- \\frac{1}{n} \\sum_n \\sum^{K}_{k=1} y_k \\frac{1}{\\hat{y_k}}) \\frac{\\partial \\hat{y}}{\\partial O_i}\\\\\n",
    "& = - \\frac{1}{n} \\sum_n (y_i \\frac{1}{\\hat{y_i}} \\hat{y_i} ( 1 - \\hat{y_i} ) + \\sum^K_{k \\not= i} y_k \\frac{1}{\\hat{y_k}}( - \\hat{y_k} \\hat{y_i}))\\\\\n",
    "& = - \\frac{1}{n} \\sum_n ( y_i - y_i \\hat{y_i} - \\sum^K_{k \\not= i} y_k \\hat{y_i})\\\\\n",
    "& = - \\frac{1}{n} \\sum_n ( y_i  - \\hat{y_i} \\sum^K_{k = 1} y_k )\n",
    "\\end{aligned}\n",
    "$$\n",
    "\n",
    "由于我们处理的多类分类任务，一个样本只对应一个标记，所以$\\sum^K_{k = 1} y_k = 1$，上式在这种问题中，即可化简为：\n",
    "\n",
    "$$\\begin{aligned}\n",
    "\\frac{\\partial \\mathrm{loss}}{\\partial O_i} &= - \\frac{1}{n} \\sum_n ( y_i  - \\hat{y_i})\\\\\n",
    "& = \\frac{1}{n} \\sum_n (\\hat{y_i} -  y_i)\n",
    "\\end{aligned}\n",
    "$$\n",
    "\n",
    "将其写成矩阵表达式：\n",
    "\n",
    "$$\\begin{aligned}\n",
    "\\frac{\\partial \\mathrm{loss}}{\\partial O} &= \\frac{1}{n} (\\hat{y} - y)\n",
    "\\end{aligned}\n",
    "$$\n",
    "\n",
    "也就是说，我们的损失函数对输出层的$K$个神经单元的偏导数为$\\mathrm{softmax}$激活值减去真值。\n",
    "\n",
    "接下来我们需要求损失函数对参数$W_2$和$b_2$的偏导数\n",
    "\n",
    "$$\n",
    "\\begin{aligned}\n",
    "\\frac{\\partial loss}{\\partial W_2} & = \\frac{\\partial \\mathrm{loss}}{\\partial \\hat{y}} \\frac{\\partial \\hat{y}}{\\partial O} \\frac{\\partial O}{\\partial W_2}\\\\\n",
    "& = \\frac{\\partial loss}{\\partial O} \\frac{\\partial O}{\\partial W_2}\\\\\n",
    "& = \\frac{1}{n} (\\hat{y} - y) \\frac{\\partial O}{\\partial W_2}\\\\\n",
    "& = \\frac{1}{n} [{H_1}^\\mathrm{T} (\\hat{y} - y)]\n",
    "\\end{aligned}\n",
    "$$\n",
    "\n",
    "$$\n",
    "\\begin{aligned}\n",
    "\\frac{\\partial loss}{\\partial b_2} & = \\frac{\\partial \\mathrm{loss}}{\\partial \\hat{y}} \\frac{\\partial \\hat{y}}{\\partial O} \\frac{\\partial O}{\\partial b_2}\\\\\n",
    "& = \\frac{\\partial loss}{\\partial O} \\frac{\\partial O}{\\partial b_2}\\\\\n",
    "& = \\frac{1}{n} (\\hat{y} - y) \\frac{\\partial O}{\\partial b_2}\\\\\n",
    "& = \\frac{1}{n} \\sum^n_{i=1} (\\hat{y_i} - y_i)\n",
    "\\end{aligned}\n",
    "$$\n",
    "\n",
    "其中，$\\frac{\\partial loss}{\\partial W_2} \\in \\mathbb{R}^{h \\times K}$，$\\frac{\\partial loss}{\\partial b_2} \\in \\mathbb{R}^{1 \\times K}$。  \n",
    "**注意，由于$b_2$是被广播成$n \\times K$的矩阵，因此实际上$b_2$对每个样本的损失都有贡献，因此对其求偏导时，要把$n$个样本对它的偏导数加和。**\n",
    "\n",
    "同理，我们可以求得$\\mathrm{loss}$对$W_1$和$b_1$的偏导数：\n",
    "\n",
    "$$\n",
    "\\begin{aligned}\n",
    "\\frac{\\partial loss}{\\partial W_1} & = \\frac{\\partial \\mathrm{loss}}{\\partial \\hat{y}} \\frac{\\partial \\hat{y}}{\\partial O} \\frac{\\partial O}{\\partial H_1} \\frac{\\partial H_1}{\\partial Z} \\frac{\\partial Z}{\\partial W_1}\\\\\n",
    "& = \\frac{\\partial loss}{\\partial O} \\frac{\\partial O}{\\partial H_1} \\frac{\\partial H_1}{\\partial Z} \\frac{\\partial Z}{\\partial W_1}\\\\\n",
    "& = \\frac{1}{n} {X}^\\mathrm{T} [(\\hat{y} - y) {W_2}^\\mathrm{T} \\frac{\\partial H_1}{\\partial Z}]\\\\\n",
    "\\end{aligned}\n",
    "$$\n",
    "\n",
    "由于我们使用的是$\\mathrm{ReLU}$激活函数，它的偏导数为：\n",
    "\n",
    "$$\\frac{\\partial \\mathrm{ReLU(x)}}{\\partial x} = \n",
    "\\begin{cases}\n",
    "0 & \\text{if } x < 0\\\\\n",
    "1 & \\text{if } x \\geq 0\n",
    "\\end{cases}\n",
    "$$\n",
    "\n",
    "所以上式为：\n",
    "\n",
    "$$\n",
    "\\frac{\\partial loss}{\\partial {W_1}_{ij}} =\n",
    "\\begin{cases}\n",
    "0 & \\text{if } {Z}_{ij} < 0\\\\\n",
    "    \\frac{1}{n} {X}^\\mathrm{T} (\\hat{y} - y) {W_2}^\\mathrm{T} & \\text{if } {Z}_{ij} \\geq 0\n",
    "\\end{cases}\n",
    "$$\n",
    "\n",
    "其中，${W_1}_{ij}$表示矩阵$W_1$第$i$行第$j$列的值，${Z}_{ij}$表示矩阵$Z$第$i$行第$j$列的值。  \n",
    "同理：\n",
    "\n",
    "$$\n",
    "\\begin{aligned}\n",
    "\\frac{\\partial loss}{\\partial b_1} & = \\frac{\\partial \\mathrm{loss}}{\\partial \\hat{y}} \\frac{\\partial \\hat{y}}{\\partial O} \\frac{\\partial O}{\\partial H_1} \\frac{\\partial H_1}{\\partial Z} \\frac{\\partial Z}{\\partial b_1}\\\\\n",
    "& = \\frac{\\partial loss}{\\partial O} \\frac{\\partial O}{\\partial H_1} \\frac{\\partial H_1}{\\partial Z} \\frac{\\partial Z}{\\partial b_1}\\\\\n",
    "& = \\frac{1}{n} (\\hat{y} - y) {W_2}^\\mathrm{T} \\frac{\\partial H_1}{\\partial Z}\\\\\n",
    "& = \\begin{cases}\n",
    "0 &\\text{if } {Z}_{ij} < 0\\\\\n",
    "\\frac{1}{n} \\sum_n (\\hat{y} - y) {W_2}^\\mathrm{T} &\\text{if } {Z}_{ij} \\geq 0\n",
    "\\end{cases}\n",
    "\\end{aligned}\n",
    "$$\n",
    "\n",
    "其中，$\\frac{\\partial loss}{\\partial W_1} \\in \\mathbb{R}^{m \\times h}$，$\\frac{\\partial loss}{\\partial b_1} \\in \\mathbb{R}^{1 \\times h}$。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 参数更新\n",
    "\n",
    "求得损失函数对四个参数的偏导数后，我们就可以使用梯度下降进行参数更新：\n",
    "$$\n",
    "W_2 := W_2 - \\alpha \\frac{\\partial \\mathrm{loss}}{\\partial W_2}\\\\\n",
    "b_2 := b_2 - \\alpha \\frac{\\partial \\mathrm{loss}}{\\partial b_2}\\\\\n",
    "W_1 := W_1 - \\alpha \\frac{\\partial \\mathrm{loss}}{\\partial W_1}\\\\\n",
    "b_1 := b_1 - \\alpha \\frac{\\partial \\mathrm{loss}}{\\partial b_1}\\\\\n",
    "$$\n",
    "其中，$\\alpha$是学习率"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "以上内容，就是一个三层感知机的前向传播与反向传播过程。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1. 导入数据"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "使用第一题的手写数字数据集"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "from time import time"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "from sklearn.datasets import load_digits\n",
    "from sklearn.model_selection import train_test_split"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "40%做测试集，60%做训练集"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "trainX, testX, trainY, testY = train_test_split(load_digits()['data'], load_digits()['target'], test_size = 0.4, random_state = 32)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[3 9 6 ... 1 0 7]\n"
     ]
    }
   ],
   "source": [
    "trainX.shape, trainY.shape, testX.shape, testY.shape\n",
    "print(trainY)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2. 数据预处理"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "使用和第一题一样的标准化处理方法"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "from sklearn.preprocessing import StandardScaler\n",
    "s = StandardScaler()\n",
    "trainX = s.fit_transform(trainX)\n",
    "testX = s.transform(testX)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "接下来还要处理输出。  \n",
    "我们的神经网络是针对每个样本，输出其分别属于$K$类的概率，我们要找最大的那个概率，对应的是哪个类。  \n",
    "我们当前的trainY和testY，每个样本都是一个类标，我们需要将其变成one_hot编码，也就是，假设当前样本的类别是3，我们需要把它变成一个长度为10的向量，其中第4个元素为1，其他元素都为0。得到的矩阵分别记为trainY_mat和testY_mat。  \n",
    "这样，模型训练完成后，会针对每个样本输出十个数，分别代表这个样本属于$0,1,...,9$的概率，那我们只要取最大的那个数的下标，就知道模型认为这个样本是哪类了。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "trainY_mat = np.zeros((len(trainY), 10))\n",
    "trainY_mat[np.arange(0, len(trainY), 1), trainY] = 1\n",
    "\n",
    "testY_mat = np.zeros((len(testY), 10))\n",
    "testY_mat[np.arange(0, len(testY), 1), testY] = 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "((1078, 10), (719, 10))"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "trainY_mat.shape, testY_mat.shape"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3. 参数初始化"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "这题和上一题的区别是，我们把参数用dict存起来"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "def initialize(h, K):\n",
    "    '''\n",
    "    参数初始化\n",
    "    \n",
    "    Parameters\n",
    "    ----------\n",
    "    h: int: 隐藏层单元个数\n",
    "    \n",
    "    K: int: 输出层单元个数\n",
    "    \n",
    "    Returns\n",
    "    ----------\n",
    "    parameters: dict，参数，键是\"W1\", \"b1\", \"W2\", \"b2\"\n",
    "    \n",
    "    '''\n",
    "    np.random.seed(32)\n",
    "    W_1 = np.random.normal(size = (trainX.shape[1], h)) * 0.01\n",
    "    b_1 = np.zeros((1, h))\n",
    "    \n",
    "    np.random.seed(32)\n",
    "    W_2 = np.random.normal(size = (h, K)) * 0.01\n",
    "    b_2 = np.zeros((1, K))\n",
    "    \n",
    "    parameters = {'W1': W_1, 'b1': b_1, 'W2': W_2, 'b2': b_2}\n",
    "    \n",
    "    return parameters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(64, 50)\n",
      "(1, 50)\n",
      "(50, 10)\n",
      "(1, 10)\n"
     ]
    }
   ],
   "source": [
    "# 测试样例\n",
    "parameterst = initialize(50, 10)\n",
    "print(parameterst['W1'].shape) # (64, 50)\n",
    "print(parameterst['b1'].shape) # (1, 50)\n",
    "print(parameterst['W2'].shape) # (50, 10)\n",
    "print(parameterst['b2'].shape) # (1, 10)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4. 前向传播"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "完成Z的计算"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "def linear_combination(X, W, b):\n",
    "    '''\n",
    "    计算Z，Z = XW + b\n",
    "    \n",
    "    Parameters\n",
    "    ----------\n",
    "    X: np.ndarray, shape = (n, m)，输入的数据\n",
    "    \n",
    "    W: np.ndarray, shape = (m, h)，权重\n",
    "    \n",
    "    b: np.ndarray, shape = (1, h)，偏置\n",
    "    \n",
    "    Returns\n",
    "    ----------\n",
    "    Z: np.ndarray, shape = (n, h)，线性组合后的值\n",
    "    \n",
    "    '''\n",
    "    \n",
    "    # Z = XW + b\n",
    "    # YOUR CODE HERE\n",
    "    Z = np.dot(X,W) + b\n",
    "    \n",
    "    return Z"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(1078, 50)\n",
      "-5.273044421225233e-19\n"
     ]
    }
   ],
   "source": [
    "# 测试样例\n",
    "parameterst = initialize(50, 10)\n",
    "Zt = linear_combination(trainX, parameterst['W1'], parameterst['b1'])\n",
    "print(Zt.shape) # (1078, 50)\n",
    "print(Zt.mean()) # -5.27304442123e-19"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$\\rm ReLU$激活函数"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "def ReLU(X):\n",
    "    '''\n",
    "    ReLU激活函数\n",
    "    \n",
    "    Parameters\n",
    "    ----------\n",
    "    X: np.ndarray，待激活的矩阵\n",
    "    \n",
    "    Returns\n",
    "    ----------\n",
    "    activations: np.ndarray, 激活后的矩阵\n",
    "    \n",
    "    '''\n",
    "    \n",
    "    # YOUR CODE HERE\n",
    "    activations = X.copy()\n",
    "    activations[activations < 0] = 0\n",
    "    return activations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.030445670920609378\n",
      "(1078, 10)\n",
      "0.0006001926584638006\n"
     ]
    }
   ],
   "source": [
    "# 测试样例\n",
    "parameterst = initialize(50, 10)\n",
    "Zt = linear_combination(trainX, parameterst['W1'], parameterst['b1'])\n",
    "Ht = ReLU(Zt)\n",
    "print(Ht.mean()) # 0.0304\n",
    "\n",
    "Ot = linear_combination(Ht, parameterst['W2'], parameterst['b2'])\n",
    "print(Ot.shape) # (1078, 10)\n",
    "print(Ot.mean()) # 0.0006"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$\\rm softmax$激活  \n",
    "\n",
    "$$\n",
    "\\mathrm{softmax}(O_i) = \\frac{\\exp{(O_i)}}{\\sum^{K}_{k=1} \\exp{(O_k)}}\n",
    "$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "def my_softmax(O):\n",
    "    '''\n",
    "    softmax激活\n",
    "    '''\n",
    "    # YOUR CODE HERE\n",
    "    softmax = np.exp(O) / np.exp(O).sum()\n",
    "    \n",
    "    return softmax"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[0.33333333 0.33333333 0.33333333]]\n",
      "[[nan nan nan]]\n",
      "[[nan nan nan]]\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "<ipython-input-15-56dcc01f6be4>:6: RuntimeWarning: invalid value encountered in true_divide\n",
      "  softmax = np.exp(O) / np.exp(O).sum()\n",
      "<ipython-input-15-56dcc01f6be4>:6: RuntimeWarning: overflow encountered in exp\n",
      "  softmax = np.exp(O) / np.exp(O).sum()\n"
     ]
    }
   ],
   "source": [
    "# 测试样例1\n",
    "print(my_softmax(np.array([[0.3, 0.3, 0.3]])))  # array([[ 0.33333333,  0.33333333,  0.33333333]])\n",
    "\n",
    "# 测试样例2\n",
    "test1 = np.array([[-1e32, -1e32, -1e32]])\n",
    "test2 = np.array([[1e32, 1e32, 1e32]])\n",
    "print(my_softmax(test1))\n",
    "print(my_softmax(test2))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "这里，其实是有数值计算上的问题的，假设，我们最后的输出有三个数，每个数都特别小，理论上来说，通过$\\rm softmax$激活后，三个值都是$\\frac{1}{3}$。但实际上就不是这样了，实际上会导致分母为0，除法就不能做了。如果每个数都特别大，会导致做指数运算的时候上溢。\n",
    "\n",
    "我们需要用其他的方法来实现$\\rm softmax$。\n",
    "\n",
    "我们将传入$\\rm softmax$的向量，每个元素减去他们中的最大值，即\n",
    "\n",
    "$$\n",
    "\\mathrm{softmax}(O_i) = \\mathrm{softmax}(O_i - \\mathrm{max(O)})\n",
    "$$\n",
    "\n",
    "这个式子是成立的，感兴趣的同学可以证明一下上面的式子。\n",
    "\n",
    "当我们做了这样的变换后，向量$O$中的最大值就变成了0，就不会上溢了，而分母中最少有一项为1，也不会出现下溢导致分母为0的问题了。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "def softmax(O):\n",
    "    '''\n",
    "    softmax激活函数\n",
    "    \n",
    "    Parameters\n",
    "    ----------\n",
    "    O: np.ndarray，待激活的矩阵\n",
    "    \n",
    "    Returns\n",
    "    ----------\n",
    "    activations: np.ndarray, 激活后的矩阵\n",
    "    \n",
    "    '''\n",
    "    \n",
    "    # YOUR CODE HEER\n",
    "    tmp = O - np.max(O)\n",
    "    activations = np.exp(tmp) / np.exp(tmp).sum(axis = 1, keepdims = True)\n",
    "    \n",
    "    return activations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(1078, 10)\n",
      "0.0006001926584638006\n",
      "0.1\n"
     ]
    }
   ],
   "source": [
    "# 测试样例\n",
    "parameterst = initialize(50, 10)\n",
    "Zt = linear_combination(trainX, parameterst['W1'], parameterst['b1'])\n",
    "Ht = ReLU(Zt)\n",
    "Ot = linear_combination(Ht, parameterst['W2'], parameterst['b2'])\n",
    "y_pred = softmax(Ot)\n",
    "\n",
    "print(y_pred.shape)  # (1078, 10)\n",
    "print(Ot.mean())     # 0.000600192658464\n",
    "print(y_pred.mean()) # 0.1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "接下来是实现损失函数，交叉熵损失函数：\n",
    "\n",
    "$$\n",
    "\\mathrm{loss} = - \\frac{1}{n} \\sum_n \\sum^{K}_{k=1} y_k \\log{(\\hat{y_k})}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "这里又会出一个问题，交叉熵损失函数中，我们需要对$\\rm softmax$的激活值取对数，也就是$\\log{\\hat{y}}$，这就要求我们的激活值全都是大于0的数，不能等于0，但是我们实现的$\\rm softmax$在有些时候确实会输出0，比如："
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[1., 0., 0.]])"
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "softmax(np.array([[1e32, 0, -1e32]]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "这就使得在计算loss的时候会出现问题，解决这个问题的方法是$\\rm log \\ softmax$。所谓$\\rm log \\ softmax$，就是将交叉熵中的对数运算与$\\rm softmax$结合起来，避开为0的情况\n",
    "\n",
    "$$\\begin{aligned}\n",
    "\\log{\\frac{\\exp{(O_i)}}{\\sum_K \\exp{(O_k)}}} &= \\log{\\frac{\\exp{(O_i - \\mathrm{max}(O))}}{\\sum_K \\exp{(O_k - \\mathrm{max}(O))}}}\\\\\n",
    "&= O_i - \\mathrm{max}(O) - \\log{\\sum_K \\exp{(O_k - \\mathrm{max}(O))}}\n",
    "\\end{aligned}\n",
    "$$\n",
    "\n",
    "这样我们再计算$\\rm loss$的时候就可以把输出层的输出直接放到$\\rm log \\ softmax$中计算，不用先激活，再取对数了。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "我们先编写`log_softmax`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "def log_softmax(x):\n",
    "    '''\n",
    "    log softmax\n",
    "    \n",
    "    Parameters\n",
    "    ----------\n",
    "    x: np.ndarray，待激活的矩阵\n",
    "    \n",
    "    Returns\n",
    "    ----------\n",
    "    log_activations: np.ndarray, 激活后取了对数的矩阵\n",
    "    \n",
    "    '''\n",
    "    # YOUR CODE HERE\n",
    "    log_activations = x - np.max(x) - np.log(np.sum(np.exp(x - np.max(x)),axis = 1, keepdims = True))\n",
    "    \n",
    "    return log_activations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(1078, 10)\n",
      "-2.3025914871652615\n"
     ]
    }
   ],
   "source": [
    "# 测试样例\n",
    "parameterst = initialize(50, 10)\n",
    "Zt = linear_combination(trainX, parameterst['W1'], parameterst['b1'])\n",
    "Ht = ReLU(Zt)\n",
    "Ot = linear_combination(Ht, parameterst['W2'], parameterst['b2'])\n",
    "t = log_softmax(Ot)\n",
    "print(t.shape)  # (1078, 10)\n",
    "print(t.mean()) # -2.30259148717"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "然后编写`cross_entropy_with_softmax`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "def cross_entropy_with_softmax(y_true, O):\n",
    "    '''\n",
    "    求解交叉熵损失函数，这里需要使用log softmax，所以参数分别是真值和未经softmax激活的输出值\n",
    "\n",
    "    Parameters\n",
    "    ----------\n",
    "    y_true: np.ndarray，shape = (n, K), 真值\n",
    "    \n",
    "    O: np.ndarray, shape = (n, K)，softmax激活前的输出层的输出值\n",
    "    \n",
    "    Returns\n",
    "    ----------\n",
    "    loss: float, 平均的交叉熵损失值\n",
    "    \n",
    "    '''\n",
    "    \n",
    "    # 平均交叉熵损失\n",
    "    # YOUR CODE HERE\n",
    "    loss = - np.sum(log_softmax(O) * y_true) / len(y_true)\n",
    "    \n",
    "    return loss"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2.302667079581974\n"
     ]
    }
   ],
   "source": [
    "# 测试样例\n",
    "parameterst = initialize(50, 10)\n",
    "Zt = linear_combination(trainX, parameterst['W1'], parameterst['b1'])\n",
    "Ht = ReLU(Zt)\n",
    "Ot = linear_combination(Ht, parameterst['W2'], parameterst['b2'])\n",
    "losst = cross_entropy_with_softmax(trainY_mat, Ot)\n",
    "print(losst.mean()) # 2.30266707958"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**正是因为$\\rm softmax$激活与交叉熵损失会有这样的问题，所以在很多深度学习框架中，交叉熵损失函数就直接带有了激活的功能，所以我们在实现前向传播计算的时候，就不要加$\\rm softmax$激活函数了。**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "def forward(X, parameters):\n",
    "    '''\n",
    "    前向传播，从输入一直到输出层softmax激活前的值\n",
    "    \n",
    "    Parameters\n",
    "    ----------\n",
    "    X: np.ndarray, shape = (n, m)，输入的数据\n",
    "    \n",
    "    parameters: dict，参数\n",
    "    \n",
    "    Returns\n",
    "    ----------\n",
    "    O: np.ndarray, shape = (n, K)，softmax激活前的输出层的输出值\n",
    "    \n",
    "    '''\n",
    "    # 输入层到隐藏层\n",
    "    # YOUR CODE HERE\n",
    "    Z = linear_combination(X, parameters['W1'], parameters['b1'])\n",
    "    \n",
    "    # 隐藏层的激活\n",
    "    # YOUR CODE HERE\n",
    "    H = ReLU(Z)\n",
    "    \n",
    "    # 隐藏层到输出层\n",
    "    # YOUR CODE HERE\n",
    "    O = linear_combination(H, parameters['W2'], parameters['b2'])\n",
    "\n",
    "    return O"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.0006001926584638006\n"
     ]
    }
   ],
   "source": [
    "# 测试样例\n",
    "parameterst = initialize(50, 10)\n",
    "Ot = forward(trainX, parameterst)\n",
    "print(Ot.mean()) # 0.000600192658464"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 5. 反向传播"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "先计算梯度"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "def compute_gradient(y_true, y_pred, H, Z, X, parameters):\n",
    "    '''\n",
    "    计算梯度\n",
    "    \n",
    "    Parameters\n",
    "    ----------\n",
    "    y_true: np.ndarray，shape = (n, K), 真值\n",
    "    \n",
    "    y_pred: np.ndarray, shape = (n, K)，softmax激活后的输出层的输出值\n",
    "    \n",
    "    H: np.ndarray, shape = (n, h)，隐藏层激活后的值\n",
    "    \n",
    "    Z: np.ndarray, shape = (n, h), 隐藏层激活前的值\n",
    "    \n",
    "    X: np.ndarray, shape = (n, m)，输入的原始数据\n",
    "    \n",
    "    parameters: dict，参数\n",
    "    \n",
    "    Returns\n",
    "    ----------\n",
    "    grads: dict, 梯度\n",
    "    \n",
    "    '''\n",
    "    \n",
    "    # 计算W2的梯度\n",
    "    # YOUR CODE HERE\n",
    "    dW2 = np.dot(H.T, (y_pred - y_true)) / len(y_pred)\n",
    "    \n",
    "    # 计算b2的梯度\n",
    "    # YOUR CODE HERE\n",
    "    db2 = np.sum(y_pred - y_true, axis = 0) / len(y_pred)\n",
    "    \n",
    "    # 计算ReLU的梯度\n",
    "    relu_grad = Z.copy()\n",
    "    relu_grad[relu_grad >= 0] = 1\n",
    "    relu_grad[relu_grad < 0] = 0\n",
    "    \n",
    "    # 计算W1的梯度\n",
    "    # YOUR CODE HERE\n",
    "    dW1 = np.dot(X.T, np.dot(y_pred - y_true, parameters['W2'].T) * relu_grad) / len(y_pred)\n",
    "    \n",
    "    # 计算b1的梯度\n",
    "    # YOUR CODE HERE\n",
    "    db1 = np.sum((np.dot(y_pred - y_true, parameters['W2'].T) * relu_grad), axis = 0) / len(y_pred)\n",
    "    \n",
    "    grads = {'dW2': dW2, 'db2': db2, 'dW1': dW1, 'db1': db1}\n",
    "    \n",
    "    return grads"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.0513287739604719\n",
      "-0.0005277166289567186\n",
      "-5.204170427930421e-18\n",
      "-1.5178830414797062e-17\n"
     ]
    }
   ],
   "source": [
    "# 测试样例\n",
    "parameterst = initialize(50, 10)\n",
    "\n",
    "Zt = linear_combination(trainX, parameterst['W1'], parameterst['b1'])\n",
    "Ht = ReLU(Zt)\n",
    "Ot = linear_combination(Ht, parameterst['W2'], parameterst['b2'])\n",
    "y_predt = softmax(Ot)\n",
    "\n",
    "gradst = compute_gradient(trainY_mat, y_predt, Ht, Zt, trainX, parameterst)\n",
    "\n",
    "print(gradst['dW1'].sum()) # 0.0429186117668\n",
    "print(gradst['db1'].sum()) # -5.05985151857e-05\n",
    "print(gradst['dW2'].sum()) # -2.16840434497e-18\n",
    "print(gradst['db2'].sum()) # -1.34441069388e-17"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "梯度下降，参数更新"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "def update(parameters, grads, learning_rate):\n",
    "    '''\n",
    "    参数更新\n",
    "    \n",
    "    Parameters\n",
    "    ----------\n",
    "    parameters: dict，参数\n",
    "    \n",
    "    grads: dict, 梯度\n",
    "    \n",
    "    learning_rate: float, 学习率\n",
    "    \n",
    "    '''\n",
    "    parameters['W2'] -= learning_rate * grads['dW2']\n",
    "    parameters['b2'] -= learning_rate * grads['db2']\n",
    "    parameters['W1'] -= learning_rate * grads['dW1']\n",
    "    parameters['b1'] -= learning_rate * grads['db1']"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "反向传播，参数更新"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.5834954544808615\n",
      "0.0\n",
      "0.18887164310031426\n",
      "0.0\n",
      "\n",
      "0.5783625770848143\n",
      "5.2771662895671786e-05\n",
      "0.18887164310031426\n",
      "1.463672932855431e-18\n"
     ]
    }
   ],
   "source": [
    "# 测试样例\n",
    "parameterst = initialize(50, 10)\n",
    "print(parameterst['W1'].sum())  # 0.583495454481\n",
    "print(parameterst['b1'].sum())  # 0.0\n",
    "print(parameterst['W2'].sum())  # 0.1888716431\n",
    "print(parameterst['b2'].sum())  # 0.0\n",
    "print()\n",
    "\n",
    "Zt = linear_combination(trainX, parameterst['W1'], parameterst['b1'])\n",
    "Ht = ReLU(Zt)\n",
    "Ot = linear_combination(Ht, parameterst['W2'], parameterst['b2'])\n",
    "y_predt = softmax(Ot)\n",
    "\n",
    "gradst = compute_gradient(trainY_mat, y_predt, Ht, Zt, trainX, parameterst)\n",
    "update(parameterst, gradst, 0.1)\n",
    "\n",
    "print(parameterst['W1'].sum())  # 0.579203593304\n",
    "print(parameterst['b1'].sum())  # 5.05985151857e-06\n",
    "print(parameterst['W2'].sum())  # 0.1888716431\n",
    "print(parameterst['b2'].sum())  # 1.24683249836e-18"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "def backward(y_true, y_pred, H, Z, X, parameters, learning_rate):\n",
    "    '''\n",
    "    计算梯度，参数更新\n",
    "    \n",
    "    Parameters\n",
    "    ----------\n",
    "    y_true: np.ndarray，shape = (n, K), 真值\n",
    "    \n",
    "    y_pred: np.ndarray, shape = (n, K)，softmax激活后的输出层的输出值\n",
    "    \n",
    "    H: np.ndarray, shape = (n, h)，隐藏层激活后的值\n",
    "    \n",
    "    Z: np.ndarray, shape = (n, h), 隐藏层激活前的值\n",
    "    \n",
    "    X: np.ndarray, shape = (n, m)，输入的原始数据\n",
    "    \n",
    "    parameters: dict，参数\n",
    "    \n",
    "    learning_rate: float, 学习率\n",
    "    \n",
    "    '''\n",
    "    # 计算梯度\n",
    "    # YOUR CODE HERE\n",
    "    grads = compute_gradient(y_true, y_pred, H, Z, X, parameters)\n",
    "    \n",
    "    # 更新参数\n",
    "    # YOUR CODE HERE\n",
    "    update(parameters, grads, learning_rate)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.5834954544808615\n",
      "0.0\n",
      "0.18887164310031426\n",
      "0.0\n",
      "\n",
      "0.5783625770848143\n",
      "5.2771662895671786e-05\n",
      "0.18887164310031426\n",
      "1.463672932855431e-18\n"
     ]
    }
   ],
   "source": [
    "# 测试样例\n",
    "parameterst = initialize(50, 10)\n",
    "print(parameterst['W1'].sum())  # 0.583495454481\n",
    "print(parameterst['b1'].sum())  # 0.0\n",
    "print(parameterst['W2'].sum())  # 0.1888716431\n",
    "print(parameterst['b2'].sum())  # 0.0\n",
    "print()\n",
    "\n",
    "Zt = linear_combination(trainX, parameterst['W1'], parameterst['b1'])\n",
    "Ht = ReLU(Zt)\n",
    "Ot = linear_combination(Ht, parameterst['W2'], parameterst['b2'])\n",
    "y_predt = softmax(Ot)\n",
    "\n",
    "backward(trainY_mat, y_predt, Ht, Zt, trainX, parameterst, 0.1)\n",
    "\n",
    "print(parameterst['W1'].sum())  # 0.579203593304\n",
    "print(parameterst['b1'].sum())  # 5.05985151857e-06\n",
    "print(parameterst['W2'].sum())  # 0.1888716431\n",
    "print(parameterst['b2'].sum())  # 1.24683249836e-18"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 6. 训练"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "def train(trainX, trainY, testX, testY, parameters, epochs, learning_rate = 0.01, verbose = False):\n",
    "    '''\n",
    "    训练\n",
    "    \n",
    "    Parameters\n",
    "    ----------\n",
    "    Parameters\n",
    "    ----------\n",
    "    trainX: np.ndarray, shape = (n, m), 训练集\n",
    "    \n",
    "    trainY: np.ndarray, shape = (n, K), 训练集标记\n",
    "    \n",
    "    testX: np.ndarray, shape = (n_test, m)，测试集\n",
    "    \n",
    "    testY: np.ndarray, shape = (n_test, K)，测试集的标记\n",
    "    \n",
    "    parameters: dict，参数\n",
    "    \n",
    "    epochs: int, 要迭代的轮数\n",
    "    \n",
    "    learning_rate: float, default 0.01，学习率\n",
    "    \n",
    "    verbose: boolean, default False，是否打印损失值\n",
    "    \n",
    "    Returns\n",
    "    ----------\n",
    "    training_loss_list: list(float)，每迭代一次之后，训练集上的损失值\n",
    "    \n",
    "    testing_loss_list: list(float)，每迭代一次之后，测试集上的损失值\n",
    "    \n",
    "    '''\n",
    "    # 存储损失值\n",
    "    training_loss_list = []\n",
    "    testing_loss_list = []\n",
    "    \n",
    "    for i in range(epochs):\n",
    "        \n",
    "        # 这里要计算出Z和H，因为后面反向传播计算梯度的时候需要这两个矩阵\n",
    "        Z = linear_combination(trainX, parameters['W1'], parameters['b1'])\n",
    "        H = ReLU(Z)\n",
    "        train_O = linear_combination(H, parameters['W2'], parameters['b2'])\n",
    "        train_y_pred = softmax(train_O)\n",
    "        training_loss = cross_entropy_with_softmax(trainY, train_O)\n",
    "        \n",
    "        test_O = forward(testX, parameters)\n",
    "        testing_loss = cross_entropy_with_softmax(testY, test_O)\n",
    "        \n",
    "        if verbose == True:\n",
    "            print('epoch %s, training loss:%s'%(i + 1, training_loss))\n",
    "            print('epoch %s, testing loss:%s'%(i + 1, testing_loss))\n",
    "            print()\n",
    "        \n",
    "        training_loss_list.append(training_loss)\n",
    "        testing_loss_list.append(testing_loss)\n",
    "        \n",
    "        backward(trainY, train_y_pred, H, Z, trainX, parameters, learning_rate)\n",
    "    return training_loss_list, testing_loss_list"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.5834954544808615\n",
      "0.0\n",
      "0.18887164310031426\n",
      "0.0\n",
      "\n",
      "0.5783625770848143\n",
      "5.2771662895671786e-05\n",
      "0.18887164310031426\n",
      "1.463672932855431e-18\n"
     ]
    }
   ],
   "source": [
    "# 测试样例\n",
    "parameterst = initialize(50, 10)\n",
    "print(parameterst['W1'].sum())  # 0.583495454481\n",
    "print(parameterst['b1'].sum())  # 0.0\n",
    "print(parameterst['W2'].sum())  # 0.1888716431\n",
    "print(parameterst['b2'].sum())  # 0.0\n",
    "print()\n",
    "\n",
    "training_loss_list, testing_loss_list = train(trainX, trainY_mat, testX, testY_mat, parameterst, 1, 0.1, False)\n",
    "\n",
    "print(parameterst['W1'].sum())  # 0.579203593304\n",
    "print(parameterst['b1'].sum())  # 5.05985151857e-06\n",
    "print(parameterst['W2'].sum())  # 0.1888716431\n",
    "print(parameterst['b2'].sum())  # 1.24683249836e-18"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 7. 绘制模型损失值变化曲线"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "def plot_loss_curve(training_loss_list, testing_loss_list):\n",
    "    '''\n",
    "    绘制损失值变化曲线\n",
    "    \n",
    "    Parameters\n",
    "    ----------\n",
    "    training_loss_list: list(float)，每迭代一次之后，训练集上的损失值\n",
    "    \n",
    "    testing_loss_list: list(float)，每迭代一次之后，测试集上的损失值\n",
    "    \n",
    "    '''\n",
    "    plt.figure(figsize = (10, 6))\n",
    "    plt.plot(training_loss_list, label = 'training loss')\n",
    "    plt.plot(testing_loss_list, label = 'testing loss')\n",
    "    plt.xlabel('epoch')\n",
    "    plt.ylabel('loss')\n",
    "    plt.legend()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 8. 预测"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "模型训练完后，我们的就可以进行预测了，需要注意的是，我们的神经网络是针对每个样本，输出其分别属于$K$类的概率，我们要找最大的那个概率，对应的是哪个类。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "def predict(X, parameters):\n",
    "    '''\n",
    "    预测，调用forward函数完成神经网络对输入X的计算，然后完成类别的划分，取每行最大的那个数的下标作为标记\n",
    "    \n",
    "    Parameters\n",
    "    ----------\n",
    "    X: np.ndarray, shape = (n, m), 训练集\n",
    "    \n",
    "    parameters: dict，参数\n",
    "    \n",
    "    Returns\n",
    "    ----------\n",
    "    prediction: np.ndarray, shape = (n, 1)，预测的标记\n",
    "    \n",
    "    '''\n",
    "    # 用forward函数得到softmax激活前的值\n",
    "    # YOUR CODE HERE\n",
    "    O = forward(X, parameters)\n",
    "    \n",
    "    # 计算softmax激活后的值\n",
    "    # YOUR CODE HERE\n",
    "    y_pred = softmax(O)\n",
    "    \n",
    "    # 取每行最大的元素对应的下标\n",
    "    # YOUR CODE HERE\n",
    "    prediction = np.argmax(y_pred,axis = 1)\n",
    "    \n",
    "    return prediction"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.15577190542420027"
      ]
     },
     "execution_count": 36,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# 测试样例\n",
    "from sklearn.metrics import accuracy_score\n",
    "\n",
    "parameterst = initialize(50, 10)\n",
    "training_loss_list, testing_loss_list = train(trainX, trainY_mat, testX, testY_mat, parameterst, 1, 0.1, False)\n",
    "\n",
    "predictiont = predict(testX, parameterst)\n",
    "accuracy_score(predictiont, testY)  # 0.15994436717663421"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 9. 训练一个三层感知机"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "隐藏层单元数设置为50，输出层单元数为10，我们设置学习率为0.03，迭代轮数为1000轮"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "training time: 3.692716598510742 s\n"
     ]
    }
   ],
   "source": [
    "start_time = time()\n",
    "\n",
    "h = 50\n",
    "K = 10\n",
    "parameters = initialize(h, K)\n",
    "training_loss_list, testing_loss_list = train(trainX, trainY_mat, testX, testY_mat, parameters, 1000, 0.03, False)\n",
    "\n",
    "end_time = time()\n",
    "print('training time: %s s'%(end_time - start_time))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "计算测试集精度"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 38,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.9457579972183588"
      ]
     },
     "execution_count": 38,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "prediction = predict(testX, parameters)\n",
    "accuracy_score(prediction, testY)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "绘制损失值变化曲线"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAmEAAAFzCAYAAAB2A95GAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nOzdd3QV1cLG4d9O7yGNGiChJJQktNB7LyoodsSOKJaLDcVr1+tVrx0LXiuKDRUbgoL03nsJCYTQAqRASO/z/XFy/Sx0EiblfdbKSjIz5/gedC1e9+zZ21iWhYiIiIhcWE52BxARERGpiVTCRERERGygEiYiIiJiA5UwERERERuohImIiIjYQCVMRERExAYudgc4W8HBwVZYWJjdMUREREROa926dWmWZYWc6FyVK2FhYWGsXbvW7hgiIiIip2WM2Xuyc7odKSIiImIDlTARERERG6iEiYiIiNigys0JExERkb8rKiriwIED5Ofn2x2lRvLw8CA0NBRXV9czfo1KmIiISDVw4MABfH19CQsLwxhjd5waxbIs0tPTOXDgAOHh4Wf8Ot2OFBERqQby8/MJCgpSAbOBMYagoKCzHoVUCRMREakmVMDscy5/9iphIiIict4yMjJ45513zum1w4YNIyMj45TXPPHEE8ydO/ec3v+vwsLCSEtLK5f3Oh8qYSIiInLeTlXCSkpKTvnaWbNmUatWrVNe88wzzzBgwIBzzlcZqYSJiIjIeZs4cSK7d++mbdu2TJgwgYULF9K3b19GjRpFdHQ0AJdeeikdOnSgdevWvPfee7+/9n8jU0lJSbRs2ZLbbruN1q1bM2jQIPLy8gC46aab+Pbbb3+//sknn6R9+/ZER0cTFxcHQGpqKgMHDqR9+/bcfvvtNG7c+LQjXq+++ipRUVFERUXx+uuvA5CTk8NFF11EmzZtiIqKYtq0ab9/xlatWhETE8ODDz543n9mejpSRESkmnl6xja2J2eW63u2qu/Hk5e0Pun5F154ga1bt7Jx40YAFi5cyOrVq9m6devvTwx+9NFHBAYGkpeXR8eOHbn88ssJCgr60/skJCTw5Zdf8v7773PVVVcxffp0Ro8e/bd/XnBwMOvXr+edd97h5Zdf5oMPPuDpp5+mX79+PPLII/z6669/Knonsm7dOj7++GNWrVqFZVl07tyZ3r17k5iYSP369Zk5cyYAx48f5+jRo3z//ffExcVhjDnt7dMzoZGwv8jOPMbGBd+yY/Vv7Nm2iiN7d5J99DClhVp3RURE5Gx06tTpT0s2TJo0iTZt2tClSxf2799PQkLC314THh5O27ZtAejQoQNJSUknfO+RI0f+7ZqlS5dyzTXXADBkyBACAgJOmW/p0qVcdtlleHt74+Pjw8iRI1myZAnR0dHMnTuXhx9+mCVLluDv74+fnx8eHh6MGTOG7777Di8vr7P94/gbjYT9xeE922m76NYTnivCmTw8yTOeFDh7UeTsRbGLFyUu3pS6+WDcvMHdF2cPX1w8fHH18sPbPwifgNp4+IWAVxB4BoDzmS/kJiIicrZONWJ1IXl7e//+88KFC5k7dy4rVqzAy8uLPn36nHBJB3d3999/dnZ2/v125Mmuc3Z2pri4GHCs13U2TnZ9REQE69atY9asWTzyyCMMGjSIJ554gtWrVzNv3jy++uor3nrrLebPn39W/7y/Ugn7i/pNo4i/eDoFOccpzM2kOC+LkrwsSguysAqyMYXZOBXl4Fycg2txLm5FWbiXpuBp5eFNPt7k4WZOPQEx13iR5+JPoVstSjwCMT7BuNQKxSu4Id7BDXHybwC+9cE7BJw0WCkiIpWfr68vWVlZJz1//PhxAgIC8PLyIi4ujpUrV5Z7hh49evD111/z8MMPM2fOHI4dO3bK63v16sVNN93ExIkTsSyL77//nqlTp5KcnExgYCCjR4/Gx8eHKVOmkJ2dTW5uLsOGDaNLly40a9bsvPOqhP2Fl48/EbFn//SFZVnkFpZwNL+YrJxscrKOk5uVQc7xdAoyUynKSqM05yjkH8W14BjuhRl4ZmXin3WIkLQdBJGBiyn903uW4EyOR12K/MNwrd0M73qROAc3g8CmEBgOTs7l9bFFRETOS1BQEN27dycqKoqhQ4dy0UUX/en8kCFDePfdd4mJiSEyMpIuXbqUe4Ynn3ySa6+9lmnTptG7d2/q1auHr6/vSa9v3749N910E506dQJgzJgxtGvXjtmzZzNhwgScnJxwdXVl8uTJZGVlMWLECPLz87Esi9dee+2885qzHbqzW2xsrLV27Vq7Y5QLy7I4lltEckYeyceyyUhNJjd1H0UZB7Eyk3HOPkxwyRHCzGHCzWH8TO7vry128iA/MBL30La4NoiBujFQNxpcPW38RCIiYpcdO3bQsmVLu2PYqqCgAGdnZ1xcXFixYgXjxo37/UGBC+FE/w6MMessy4o90fUaCbORMYZAbzcCvd2IauAPNAA6/umaYzmFJKZl82tKNocPHyTr4E6KU+NpUJBIqyN7aZX6LbU2fgJAqXGhuF473MK7QaNu0KizYw6aiIhIDbBv3z6uuuoqSktLcXNz4/3337c70imphFVyAd5udPAOpEPjQKAR0BWAlMx8tiVnMvVABnuT4inav4GWxTvpeCCOmOS3cV32BpZxwmoQi1PEYIgYDHWiQFtaiIhINdW8eXM2bNhgd4wzphJWRdX286C2nwd9W9QGIigpvYjtyZms2pPOBwkHyd6zmlhrCwP2b6L1gWdh/rNYvvUxUSMh+kqo10aFTERExEYqYdWEs5MhOtSf6FB/xvRsQl5hN5btSmPqjiNs2B5HTP4aLspaR4+V7+Ky4i2s4AhMm2uh/Q3gHWx3fBERkRpHJaya8nRzZkCrOgxoVYeSy6JZtqsv360/wMRtCfQtXcl1GSuJmvc01sLnMa0vg05jIfSE8wZFRESkAqiE1QDOToZeESH0igghKz+KHzd2Zvyy4Zi0nYz1nM+l237GbfM0aNwDek+A8N66VSkiIlLBtBJoDePr4croLo357b7ePH7zSGY0uJ92OZN42dxEzqGd8OkI+GgwJC2zO6qIiFQhGRkZvPPOO+f8+tdff53c3P9fimnYsGHlsj9jUlISUVFR5/0+FUElrIZycjL0jghh6q2dmXpnfzaFjqJ95kv8x2UseWl7YcowmDYa0nfbHVVERKqA8i5hs2bNolatWuURrdJSCRPaNwpg6q2d+WRsLxb6Dqfdsef52vdGShPmwdudYcG/objA7pgiIlKJTZw4kd27d9O2bVsmTJgAwEsvvUTHjh2JiYnhySefBCAnJ4eLLrqINm3aEBUVxbRp05g0aRLJycn07duXvn37AhAWFkZaWhpJSUm0bNmS2267jdatWzNo0KDf95Ncs2YNMTExdO3alQkTJpx2xCs/P5+bb76Z6Oho2rVrx4IFCwDYtm0bnTp1om3btsTExJCQkHDCnOVNc8Lkd12aBDHjnh58sXofz8324fXCzkxpMIOIRS/Cth9gxFvQsJPdMUVE5HR+mQiHt5Tve9aNhqEvnPT0Cy+8wNatW39foX7OnDkkJCSwevVqLMti+PDhLF68mNTUVOrXr8/MmTMBx56S/v7+vPrqqyxYsIDg4L8/sZ+QkMCXX37J+++/z1VXXcX06dMZPXo0N998M++99x7dunVj4sSJp/0Ib7/9NgBbtmwhLi6OQYMGER8fz7vvvsv48eO57rrrKCwspKSkhFmzZv0tZ3nTSJj8ibOT4foujZn/QG/atmrBoH3X84zf0xTlZ8GHg2DuU1BSZHdMERGp5ObMmcOcOXNo164d7du3Jy4ujoSEBKKjo5k7dy4PP/wwS5Yswd/f/7TvFR4eTtu2bQHo0KEDSUlJZGRkkJWVRbdu3QAYNWrUad9n6dKlXH/99QC0aNGCxo0bEx8fT9euXfn3v//Niy++yN69e/H09DynnGdLI2FyQkE+7rw9qj0zNh/i8R9c+ank30wP/5nGS19zTNq/4kOo1cjumCIiciKnGLG6UCzL4pFHHuH222//27l169Yxa9YsHnnkEQYNGsQTTzxxyvdyd3f//WdnZ2fy8vI4l72vT/aaUaNG0blzZ2bOnMngwYP54IMP6Nev31nnPFsaCZOTMsYwvE19Zt/bi7B6dei9cyTfhj+NlbID3u0B8XPsjigiIpWEr68vWVlZv/8+ePBgPvroI7KzswE4ePAgKSkpJCcn4+XlxejRo3nwwQdZv379CV9/OgEBAfj6+rJy5UoAvvrqq9O+plevXnz++ecAxMfHs2/fPiIjI0lMTKRJkyb84x//YPjw4WzevPmkOcuTRsLktOr6e/Dl2C48PyuOB5fBwtDXed3pNVy+vBoG/Qu63Kl1xUREarigoCC6d+9OVFQUQ4cO5aWXXmLHjh107erY89jHx4fPPvuMXbt2MWHCBJycnHB1dWXy5MkAjB07lqFDh1KvXr3fJ8yfzocffshtt92Gt7c3ffr0Oe0twzvvvJM77riD6OhoXFxcmDJlCu7u7kybNo3PPvsMV1dX6tatyxNPPMGaNWtOmLM8mXMZzrNTbGystXbtWrtj1Fg/bUrmwW820ayW4ds6U/Da/Qu0vxGGvQwubnbHExGpsXbs2EHLli3tjnFBZWdn4+PjAzgeDDh06BBvvPGGbXlO9O/AGLPOsqwTbkmjkTA5K8Pb1KeevwdjPllL76RbmNmmKbXXvwWZB+GqqeDmZXdEERGpIWbOnMnzzz9PcXExjRs3ZsqUKXZHOiuaEyZnrWNYINPHdcXN1YW+G3qS1O3fsGsefH4F5GfaHU9ERGqIq6++mo0bN7J161ZmzpxJSEiI3ZHOikqYnJNmtX357s5u1PHz4KKlTdnd6w3Yvwo+HQ65R+2OJyIiUumphMk5q+PnmLBfx8+D4QvrkND3XTiyHaZeCvnlv6idiIicWlWb512dnMufvUqYnJc/FrFLf/MlacB/HUXs86ugMMfueCIiNYaHhwfp6ekqYjawLIv09HQ8PDzO6nV6OlLKxeHj+Vw+eTkFxaX8OjCd4F/vgPDecO1X4Hp2/1GKiMjZKyoq4sCBA+Tn59sdpUby8PAgNDQUV1fXPx0/1dORKmFSbnalZHH55BUEersxo+defH75B7QaAVdMAScNuoqISM1zqhKmvxml3DSr7cuHN8aSnJHH6LXNKOr/NGz/EeY/Y3c0ERGRSkclTMpVbFggb1zTjo37M3g4uTdW+5tg6Wuwfqrd0URERCoVlTApd0Oi6nLfgAi+25DMlFp3QZO+8PO9kLjI7mgiIiKVhkqYVIh7+jVjUKs6/OvXXayKfQ2CmsE3N0HGPrujiYiIVAoqYVIhnJwMr17dlibB3tzxbQJHhn0IpcXw9Q1QpCd3REREVMKkwvi4u/DeDbEUFpdyz+wsSkZMhuQN8MtDdkcTERGxnUqYVKjwYG+euyya1UlHeeNgBPR8ANZ/Aus/tTuaiIiIrVTCpMJd2q4BV3QI5c35CSxvdDs06QOzHoLUnXZHExERsY1KmFwQTw9vTXiwN/d9s4WMwW+Cmxd8e6vmh4mISI2lEiYXhLe7C5OuaUd6diFPzE+HSyfDkS0w72m7o4mIiNhCJUwumKgG/vyjf3N+2pTMLwUx0Ol2WPkOJPxmdzQREZELTiVMLqhxfZoS3cCfR3/YSlq3R6F2a/j+DshOsTuaiIjIBaUSJheUq7MTr1zVhuz8Yh7/eRfW5R9AQRbMvB+q2GbyIiIi56PCSpgxpqExZoExZocxZpsxZvwJrjHGmEnGmF3GmM3GmPYVlUcqj4g6vtw3MIJfth7mp0P+0PefsGMGbJ1udzQREZELpiJHwoqBByzLagl0Ae4yxrT6yzVDgeZlX2OByRWYRyqRsb2a0LZhLZ6esZ1jbW6HBrEw60HIOmJ3NBERkQuiwkqYZVmHLMtaX/ZzFrADaPCXy0YAn1oOK4Faxph6FZVJKg9nJ8MLl0eTmVfE87PjHU9LFubqtqSIiNQYF2ROmDEmDGgHrPrLqQbA/j/8foC/FzWMMWONMWuNMWtTU1MrKqZcYC3q+jGmZxO+XnuAlVlB0O8xiPsZtnxrdzQREZEKV+ElzBjjA0wH7rUsK/Ovp0/wkr8Ng1iW9Z5lWbGWZcWGhIRUREyxyfj+zQkN8OTR77dQ0PEOCO3k2FsyJ93uaCIiIhWqQkuYMcYVRwH73LKs705wyQGg4R9+DwWSKzKTVC6ebs48e2kUu1Nz+O+SvTB8EhRkwpzH7I4mIiJSoSry6UgDfAjssCzr1ZNc9hNwQ9lTkl2A45ZlHaqoTFI59Y2szcUx9XhrwS72ODWC7uNh0xeQuMjuaCIiIhWmIkfCugPXA/2MMRvLvoYZY+4wxtxRds0sIBHYBbwP3FmBeaQSe+KSVrg5O/Hsz9uh1wQICIef79PekiIiUm25VNQbW5a1lBPP+frjNRZwV0VlkKqjtq8H4/s357lZO5i/O5N+F78GUy+FJa9Av0ftjiciIlLutGK+VBo3dgujSYg3z8zYTkHjXhBzNSx9DVLi7I4mIiJS7lTCpNJwc3HiqUtak5Sey0dLk2Dwv8Hdx3FbUmuHiYhINaMSJpVKr4gQBraqw5vzEzhS4gMDnoZ9y2HLN3ZHExERKVcqYVLpPH5RK4pLLZ6ftQPaXQ/12zuWrMj/6zJzIiIiVZdKmFQ6jYK8uL1XE37YmMzafRkw7GXIPgKL/2N3NBERkXKjEiaV0rg+Tanj586/Zu7AatDeMSK2cjKk7rQ7moiISLlQCZNKycvNhQcGRbJxfwY/bz4EA54CN2/HlkaapC8iItWASphUWpe3D6VlPT9e/DWOfLcA6PsYJC6EHT/ZHU1EROS8qYRJpeXsZHh0WEsOHMvj0xVJEHsL1ImGX/8Jhbl2xxMRETkvKmFSqfVoHkzfyBDenL+Lo/mlMOwlyDwAy163O5qIiMh5UQmTSu+fw1qSU1DMpHkJ0LgrtB4JyybB8YN2RxMRETlnKmFS6TWv48s1nRrx2cq9JKZmw8CnwSqFec/YHU1EROScqYRJlXDfgAjcXZx44Zc4qNUIut4Fm7+Cg+vsjiYiInJOVMKkSgjxdWdcn6bM2X6EVYnp0PN+8A5xTNLXkhUiIlIFqYRJlXFrjybU9fPgxV/jsNx8oN9jsH8lbP/B7mgiIiJnTSVMqgxPN2fGD2jO+n0Z/Lb9iGMV/TpR8NuTUJRvdzwREZGzohImVcqVHUJpEuzNS7N3UoITDH4OMvbCqnftjiYiInJWVMKkSnFxduLBwZEkpGTz3foD0KQPRAyFxS9Ddqrd8URERM6YSphUOUOj6hIT6s/rcxPILyqBQf+C4jxY8Jzd0URERM6YSphUOcYYHh7SgoMZeXy2ci8EN4OOt8H6TyAlzu54IiIiZ0QlTKqk7s2C6dk8mLcX7CIrvwh6TQA3H5j3tN3RREREzohKmFRZDw1uwbHcIt5fnAjeQdDjXtg5C/ausDuaiIjIaamESZUVHerPRTH1+GDpHlKzCqDzOPCtD789rgVcRUSk0lMJkyrtwUGRFBSX8tb8BHDzgr6PwIE1sGOG3dFEREROSSVMqrTwYG+u7tiQL1bvY196LrQZBSEtYO5TUFJkdzwREZGTUgmTKm98/+Y4OxlemxsPzi4w4Ck4utvxtKSIiEglpRImVV4dPw9u7BrGDxsPknAkCyKGQKNusPBFKMi2O56IiMgJqYRJtXB776Z4uTrz+twEMAYGPgM5KbDiLbujiYiInJBKmFQLgd5u3NojnJlbDrH14HFo2BFaDodlkyA7xe54IiIif6MSJtXGrT2b4O/pymu/xTsO9H8SivNh0Yv2BhMRETkBlTCpNvw9XRnbqwnz4lJYv++YYzujDjfBuimQvtvueCIiIn+iEibVyk3dwgjyduPVOWWjYX0mgrM7zHvG3mAiIiJ/oRIm1Yq3uwvj+jRl6a40VuxOB5/a0PUu2P4DJG+0O56IiMjvVMKk2hndpTF1/Tx4Zc5OLMuCbneDZwDMf9buaCIiIr9TCZNqx8PVmbv7NWPt3mMsik8FD3/ocR/smgtJy+yOJyIiAqiESTV1VWxDQgM8eWVOvGM0rONt4FPXMTdMm3uLiEgloBIm1ZKbixPj+zdny8HjzNl+xLG5d++HYP9KSJhjdzwRERGVMKm+LmvXgCYh3rw6J56SUgvaXQ8BYTDvWSgttTueiIjUcCphUm25ODtx34AIdh7J4ufNyeDiBn0fhSNbYNt3dscTEZEaTiVMqrWLouvRoq4vr89NoLikFKKugNqtYcFzUFJkdzwREanBVMKkWnNyMtw/MII9aTl8t+EgODlBv8fgaCJs/NzueCIiUoOphEm1N7BVHWJC/XljbgKFxaUQORRCO8LCF6Eoz+54IiJSQ6mESbVnjOGBQZEczMhj2pp9YAz0fwKykmHNh3bHExGRGkolTGqEXs2D6RQWyFsLdpFfVALhvaBJX1jyCuRn2h1PRERqIJUwqRGMMdw/KIIjmQV8tnKv42D/xyHvKKx4295wIiJSI6mESY3RpUkQPZoFM3nhbnIKiqFBB2h5Cax4C3LS7Y4nIiI1jEqY1Cj3D4ogPaeQKcuTHAf6PgZFubD0VVtziYhIzaMSJjVK+0YB9G9Rm/cWJ5KZXwS1W0DMNbDmA8g8ZHc8ERGpQVTCpMa5b2AEx/OK+HDJHseB3g9BabFjkr6IiMgFohImNU5UA3+GRtXlw6V7OJZTCIHh0G40rJsCGfvsjiciIjWESpjUSPcNjCCnsJj/Lk50HOg1wbF+2OKX7A0mIiI1hkqY1EgRdXwZ0aY+U5bvISUrH/xDocPNsOFzSN9tdzwREakBVMKkxho/IIKiEovJC8tKV8/7wdkVFv3H3mAiIlIjqIRJjRUe7M0V7UP5fNU+Dh3PA9+60Ok22DwNUnfaHU9ERKo5lTCp0e7p3wzLsnhr/i7Hge73gqsXLHze3mAiIlLtqYRJjRYa4MU1HRsxbc1+9h/NBe9g6DIOtn0Ph7fYHU9ERKoxlTCp8e7u1wxnJ8Mb8xIcB7rdDe7+sECjYSIiUnFUwqTGq+PnwfVdGvPd+gPsTs0GzwBHEds5Ew6uszueiIhUUyphIsAdfZri4erMG3PLRsM63+EoYwv+bW8wERGptiqshBljPjLGpBhjtp7kfB9jzHFjzMayrycqKovI6QT7uHNTtzBmbE4m7nAmePg5Junvmgv7VtodT0REqqGKHAmbAgw5zTVLLMtqW/b1TAVmETmtsb2a4OPmwmu/xTsOdLoNvGvD/H/ZG0xERKqlCithlmUtBo5W1PuLlLdaXm6M6dmE2duOsOXAcXDzdizgmrQEEhfZHU9ERKoZu+eEdTXGbDLG/GKMaW1zFhFu6RFGLS9XXv2tbLHWDjeDb31Y8BxYlr3hRESkWrGzhK0HGluW1QZ4E/jhZBcaY8YaY9YaY9ampqZesIBS8/h6uHJ7r6Ys2JnKur1HwdUDej0I+1c55oeJiIiUE9tKmGVZmZZlZZf9PAtwNcYEn+Ta9yzLirUsKzYkJOSC5pSa58ZujQn2ceOVOWVzw9pdD7UaOeaGaTRMRETKiW0lzBhT1xhjyn7uVJYl3a48Iv/j5ebCnX2asXx3Ost3pYGLG/R+GA5thLiZdscTEZFqoiKXqPgSWAFEGmMOGGNuNcbcYYy5o+ySK4CtxphNwCTgGsvSMINUDqM6N6Kunwev/BaPZVkQcw0ENnXMDSsttTueiIhUAxX5dOS1lmXVsyzL1bKsUMuyPrQs613Lst4tO/+WZVmtLctqY1lWF8uylldUFpGz5eHqzD39m7Fu7zEWxaeCswv0eQRStsP27+2OJyIi1YDdT0eKVFpXdmhIaIAnr8wpGw2LGgkhLR17SpYU2x1PRESqOJUwkZNwc3FifP/mbDl4nDnbj4CTM/R9BNITYMs3dscTEZEqTiVM5BQua9eAJsHevDonntJSC1pcAnVjYNELUFJkdzwREanCVMJETsHF2Yl7B0aw80gWP285BE5O0PdROJYEGz6zO56IiFRhKmEip3FxdD0i6/jy+tx4iktKIWIwhHaExS9BUb7d8UREpIpSCRM5DScnw30DI0hMzeGHjclgDPR7DDIPwrqP7Y4nIiJVlEqYyBkY3LoOUQ38eGNePIXFpdCkD4T1hCWvQGGO3fFERKQKUgkTOQPGGB4YFMn+o3l8s26/42D/JyAnFVa9a284ERGpklTCRM5Qn4gQOjQO4M15u8gvKoGGnaD5YFj2BuRl2B1PRESqGJUwkTNkjOGBgREczszny9X7HAf7PQb5x2HFW/aGExGRKkclTOQsdGsWTNcmQby9YDe5hcVQLwZaXQorJ0NOmt3xRESkClEJEzlLDwyKIC27gE9X7HUc6PsoFOXC0tfsDSYiIlWKSpjIWYoNC6RPZAjvLtpNVn4RhERAzDWw+n3ITLY7noiIVBEqYSLn4IGBkWTkFvHxsiTHgd4PgVXiWMBVRETkDKiEiZyD6FB/BrWqw/uLEzmWUwiB4dD+Blj/KRzdY3c8ERGpAlTCRM7Rg4MjyS4sZvKi3Y4DvSaAkwssetHeYCIiUiWohImco4g6voxsF8qU5UkcOp4HfvWh4xjYPA1Sd9odT0REKjmVMJHzcN/A5mDB678lOA70uA9cvWDBc/YGExGRSk8lTOQ8hAZ4MbpLY75Zt59dKdngHQxdxsH2H+HQJrvjiYhIJaYSJnKe7urbFE9XZ16ZU3YLsuvd4OEP8zUaJiIiJ6cSJnKegnzcua1XE37ZepiN+zPAsxZ0Hw8Js2H/arvjiYhIJaUSJlIOxvRsQpC3Gy/+EodlWdD5DvAOgXnPgGXZHU9ERCohlTCRcuDj7sLd/ZqxIjGdJQlp4OYNPR+ApCWwe77d8UREpBJSCRMpJ6M6NyI0wJP/zI6jtNSC2FugViOY+ySUltodT0REKhmVMJFy4u7izP0DI9h6MJNZWw+Bizv0exwOb4Gt39odT0REKhmVMJFyNKJtAyLr+PLy7J0UlZRC1BVQNwbmPwvFBXbHExGRSuSMSpgxZrwxxs84fGiMWW+MGVTR4USqGmcnw4TBkSSl5/L12v3g5AQDn4aMfbDmA7vjiVupeAYAACAASURBVIhIJXKmI2G3WJaVCQwCQoCbgRcqLJVIFda/ZW1iGwfwxtwE8gpLoGk/aNIHFr8E+cftjiciIpXEmZYwU/Z9GPCxZVmb/nBMRP7AGMPDQ1uQklXAR8v2OA4OeBryjsHS1+0NJyIilcaZlrB1xpg5OErYbGOML6DHvUROomNYIANa1mbywt2kZxdA/bYQfSWsnAyZyXbHExGRSuBMS9itwESgo2VZuYArjluSInISE4e2JK+ohDfmlW3u3e8xKC2Ghc/bG0xERCqFMy1hXYGdlmVlGGNGA48BmtwicgrNavswqlMjPl+1z7G5d0AYdBwDGz6DlDi744mIiM3OtIRNBnKNMW2Ah4C9wKcVlkqkmhg/oDmers688EtZ6eo1AVy9Yd7T9gYTERHbnWkJK7YsywJGAG9YlvUG4FtxsUSqh2Afd8b1acrcHUdYsTsdvIOgx3jYOQuSltkdT0REbHSmJSzLGPMIcD0w0xjjjGNemIicxq09wqnv78Fzs7Y7tjPqchf4NYDZj2g7IxGRGuxMS9jVQAGO9cIOAw2AlyoslUg14uHqzIQhkWw9mMmPmw6Cm5djyYpDm2DTl3bHExERm5xRCSsrXp8D/saYi4F8y7I0J0zkDI1o04CoBn689OtO8otKIPoKaBDrmBtWkG13PBERscGZblt0FbAauBK4ClhljLmiIoOJVCdOToZ/DmtJ8vF8xwKuxsCQFyD7CCzTAq4iIjXRmd6OfBTHGmE3WpZ1A9AJeLziYolUP92aBjOgZW3eWbCbtOwCaNjRsYDr8jcde0uKiEiNcqYlzMmyrJQ//J5+Fq8VkTL/W8D1td/iHQf6P+n4Pvcp2zKJiIg9zrRI/WqMmW2MuckYcxMwE5hVcbFEqqdmtX0Y3bkRX67ex/bkTKjVELr9A7ZOh32r7I4nIiIX0JlOzJ8AvAfEAG2A9yzLergig4lUV/cNjMDf05WnZmzDsizoPh586mrJChGRGuaMbylaljXdsqz7Lcu6z7Ks7ysylEh1VsvLjQcHR7J6z1FmbjkE7j4w4Ek4uA42f2V3PBERuUBOWcKMMVnGmMwTfGUZYzIvVEiR6uaajo1oVc+Pf8/cQW5hMcRc41iy4rcnIC/D7ngiInIBnLKEWZbla1mW3wm+fC3L8rtQIUWqG2cnw1PDW5N8PJ93F+4GJye46BXISYOFz9sdT0RELgA94Shik07hgQxvU593Fyey/2gu1G8LsbfA6vfg8Ba744mISAVTCROx0SPDWuBsDM/N3OE40O8x8AyAWRPAsuwNJyIiFUolTMRG9fw9uatvU37ddphlu9LAKxAGPAX7VsDmaXbHExGRCqQSJmKzMT2b0CjQiyd/2kZhcSm0He2YpD/nccg/bnc8ERGpICphIjbzcHXm6eGt2ZWSzftLEssm6b8MOamwQJP0RUSqK5UwkUqgb4vaDI2qy6R5CexLz4X67SD2Zlj9X03SFxGpplTCRCqJJy5phYuT4YmftjpW0u/3OHgGwozxUFpidzwRESlnKmEilUQ9f0/uHxTJwp2p/Lr1sGOS/pAXHCvpr/nA7ngiIlLOVMJEKpEbuzamdX0/npqxjeyCYoi+ApoNgHnPwPEDdscTEZFypBImUom4ODvx3GXRpGQV8OqceDDGsZK+VQozH9DaYSIi1YhKmEgl07ZhLUZ3bsyU5XvYcuA4BIRB339C/K+w/Ue744mISDlRCROphB4cHEmwjzsPTd9MUUkpdB4H9drALw9pg28RkWpCJUykEvL3dOW5y6LZcSjTscG3swtcMsmxdtjcJ+2OJyIi5UAlTKSSGtiqDpe0qc+k+QnEH8lybPDd9S5YNwV2L7A7noiInCeVMJFK7KlLWuHr4cpD326mpNSCvo9CUHP46R7Iz7Q7noiInIcKK2HGmI+MMSnGmK0nOW+MMZOMMbuMMZuNMe0rKotIVRXk486Tl7Ri4/4MPl62B1w94dLJkHkQfnvc7ngiInIeKnIkbAow5BTnhwLNy77GApMrMItIlTW8TX0GtKzNS7N3kpSWAw07Qrd7HLcld82zO56IiJyjCithlmUtBo6e4pIRwKeWw0qgljGmXkXlEamqjDH869Jo3JydeGj6ZkpLLejzTwiOLLstedzuiCIicg7snBPWANj/h98PlB0Tkb+o6+/B45e0YvWeo3y4dA+4ejhuS2YdgtmP2h1PRETOgZ0lzJzg2AmXAzfGjDXGrDXGrE1NTa3gWCKV05UdQhnUqg4vzd5J3OFMCO0A3e+FDVMhfrbd8URE5CzZWcIOAA3/8HsokHyiCy3Les+yrFjLsmJDQkIuSDiRysYYw/Mjo/HzdOHerzZSUFwCfSZC7dbw412Qrf9BERGpSuwsYT8BN5Q9JdkFOG5Z1iEb84hUekE+7rx4eQxxh7N49bd4cHGHyz9wLFfx413aW1JEpAqpyCUqvgRWAJHGmAPGmFuNMXcYY+4ou2QWkAjsAt4H7qyoLCLVSf+Wdbi2UyPeW5zIqsR0qNMKBj0LCbNhzQd2xxMRkTNkrCr2f86xsbHW2rVr7Y4hYqucgmKGTVpCcYnFr/f2xNfdBT6/ApKWwthFULuF3RFFRAQwxqyzLCv2ROe0Yr5IFeTt7sKrV7Xl0PE8Hvthq+OJlhHvgJsPTB8DxQV2RxQRkdNQCROpojo0DuC+ARH8uDGZb9YeAN86MOJtOLIF5j1jdzwRETkNlTCRKuzOvs3o1jSIJ37a6tjkO3IIdBwDK96C+Dl2xxMRkVNQCROpwpydDK9f3RYfdxfu/mI9eYUlMOg5qBsN34+FjP2nfxMREbGFSphIFVfbz4PXrm5LQko2T8/Y5lhN/8pPoKQYvr0ZigvtjigiIiegEiZSDfRsHsK43k35as1+ftx4EIKawoi34MAamPuk3fFEROQEVMJEqon7B0YQ2ziAf363hV0p2dD6Uuh0O6x8B7b/ZHc8ERH5C5UwkWrCxdmJN0e1w8PVmbFT15KVX+RYxLV+e8dq+kcT7Y4oIiJ/oBImUo3U8/fk7evaszc9lwe+3kSpkxtcOQWME0y7Hgpz7I4oIiJlVMJEqpkuTYL457CWzNl+hMmLdkNAY7jiQ0jZDj/cqf0lRUQqCZUwkWrolu5hjGhbn5fn7GThzhRoNgAGPAXbf4Alr9gdT0REUAkTqZaMMTw/MprIOr6M/2oj+9Jzods/IOoKmP8viJ9td0QRkRpPJUykmvJyc+G/13fAsizGfLqGzIJiGP6mYyHX6WMgNd7uiCIiNZpKmEg11jjIm8mjO5CYmsPdX2yg2NkDrvkCnN3gq1GQd8zuiCIiNZZKmEg1171ZMM9eGsXi+FSe/Xk71GoIV30Kx5IcT0xqRX0REVuohInUANd2asRtPcP5ZMVePlmeBGHdYcTbkLQEZvxDT0yKiNjAxe4AInJhTBzakj1puTw9YxuNgrzo2+ZqyNgLC56DgDDoM9HuiCIiNYpGwkRqCGcnwxvXtKVFXT/u+WID25KPQ68J0PY6WPg8bPzC7ogiIjWKSphIDeLt7sKHN8Xi5+HCjR+tYd/RPLj4dQjvDT/dA4mL7I4oIlJjqISJ1DD1/D359NZOFJeWcv1Hq0jNsxwT9YOaw1fXwcH1dkcUEakRVMJEaqBmtX358MaOHMnM5+Ypq8l28oHR08EzAD67HFJ32h1RRKTaUwkTqaE6NA5g8nUd2HEoi9unrqXAuy7c8AM4ucDUyyBjn90RRUSqNZUwkRqsb4va/OfyGJbtSue+aRsprhUO138Hhdnw6aWQnWJ3RBGRakslTKSGu7xDKI9d1JJZWw7z4DebKKkdBaO+gaxDMHWkVtUXEakgKmEiwpieTZgwOJIfNiYzcfpmSkM7wdWfQdpOx4iYipiISLlTCRMRAO7q24zx/ZvzzboDPP7jVqym/RxFLGU7fDoCco/aHVFEpFpRCROR3907oDnj+jTl81X7eHrGdqzmg+DqzyFlh4qYiEg5UwkTkd8ZY3hocCS39ghnyvIknvppG6XNBsI1XziWrfh0uIqYiEg5UQkTkT8xxvDYRS1/3/B74nebKWk6oKyIxcMnl0DWEbtjiohUeSphIvI3xhj+Oawl/+jfnK/XHuDeaRspatIPrv0SjibCR4Ph6B67Y4qIVGkqYSJyQsYY7h8YwcShLZixKZk7P19PQVgfuOEnyM9wFLHDW+2OKSJSZamEicgp3dG7Kc+MaM1v248w5pO1ZNduBzf/CsYZpgyDvSvsjigiUiWphInIad3QNYyXr2zD8t3pXP3fFaR4hsGts8E7BKZeCnGz7I4oIlLlqISJyBm5okMoH9wYy560HEa+s5zdRYFwy2yo3Qq+GgXL3wLLsjumiEiVoRImImesb2RtvhrbhfyiEi6fvJx1ac5w00xoeQnMeRR+vg9KiuyOKSJSJaiEichZiQmtxfRx3ajl6cqo91fyy87jcOUn0ON+WPcxfH4l5GXYHVNEpNJTCRORs9Y4yJvp47rRqr4f4z5fz6QFu7H6PwEj3oakpfDhQEhLsDumiEilphImIuckyMedL2/rwsh2DXj1t3ju/nIDea2vhRt+gNx0eK8v7Jhhd0wRkUpLJUxEzpmHqzOvXNWGiUNbMGvLIa7873IOBXSA2xdDSARMGw1zn4KSYrujiohUOiphInJejDHc0bspH9wQy57UHC55cxkr0z3h5l+gw02w9DX4bCTkpNkdVUSkUlEJE5Fy0b9lHb6/qzt+Hi5c98Eq3l12AOvi12H4W7BvJUzuDomL7I4pIlJpqISJSLmJqOPLj3d3Z3DrOrzwSxxjp67jeMtrYMxccPeFT0fA3Ke1jIWICCphIlLOfD1ceXtUe564uBUL4lIY/tZStpY2htsXQbvRsPRV+GgIHEuyO6qIiK1UwkSk3BljuKVHONNu70JBUSmXvbOM91ceofSSN+GKjx3LV0zuAeunapV9EamxVMJEpMJ0aBzIL+N70q9FbZ6btYMbP17NkUbD4I4lUK8N/HS3Y3HX4wftjioicsGphIlIhQrwduPd0R14fmQ0a5OOMeT1xcxJdocbZ8DQ/zgWd32nK2z4XKNiIlKjqISJSIUzxnBtp0b8/I8eNAjwZOzUdUz8fitZbW6BccugTiv48U744mrI2G93XBGRC0IlTEQumKYhPnw3rjt39G7K12v3M/i1xSxK94ObZsHg52HPYni7Myx/U09Qiki1pxImIheUm4sTE4e2YPq4bni5u3DjR6t5+LutZLa7De5aBWE9YM5j8F4f2L/a7rgiIhVGJUxEbNGuUQA/39ODO/s05Zt1+xn06mLmHfaAUdPg6s8g75hjI/AZ4yH3qN1xRUTKnUqYiNjGw9WZh4a04Ps7u+Pn6cKtn6zl9s/WkVxvANy1Grre7VjGYlJbWPEOFBfaHVlEpNwYq4o9jRQbG2utXbvW7hgiUs4Ki0v5YGkik+Yl4GQM9w5ozs3dw3FN3wmzH4Xd8yCwKQx+DiKGgDF2RxYROS1jzDrLsmJPdE4jYSJSKbi5OHFnn2b8dl9vujUN4t+z4rh40lLW5taB67+D674FJ2f48hrH9keHt9odWUTkvKiEiUil0jDQiw9u7Mh713cgK7+IK95dwb1fbSA5pAeMWw5DX4LDm+HdHjB9DKTvtjuyiMg50e1IEam0cguLeXvBLt5fsgcnA2N7NuH23k3xLs2CZZNg1btQXADtr4deD4F/A7sji4j8yaluR6qEiUild+BYLi/+upMZm5Kp7evOhMGRXN4+FKecFFjyCqz9CIwTdBwDPe8H72C7I4uIACphIlJNrNt7jGd/3s7G/RlENfDjkaEt6d4sGI7thUX/gU1fgIsnxN7seLLSr57dkUWkhlMJE5Fqo7TUYsbmZF78JY7k4/l0axrEg4Mjad8oAFLjYcnLsKVsEn+70dD9XghobHdsEamhVMJEpNrJLyrhy9X7eGv+LtJzChnQsg4PDo6gRV0/OLoHlr0OG7+A0hKIuQp63A8hEXbHFpEaRiVMRKqtnIJiPl62h/8uTiS7oJgRberzj/7NaRLiA5nJjn0o134MxfnQ4iLoehc06qp1xkTkgrCthBljhgBvAM7AB5ZlvfCX832AH4E9ZYe+syzrmVO9p0qYiJxIRm4h7y1O5ONlSRQUl3BxTH3u6tuMyLq+kJMGKyfD2g8d2yHVa+uYM9b6UnB2tTu6iFRjtpQwY4wzEA8MBA4Aa4BrLcva/odr+gAPWpZ18Zm+r0qYiJxKWnYBHyzZw9QVSeQUljC4dR3u6decqAb+UJgLm76Ele9A+i7wrQ+dboMON4FXoN3RRaQasquEdQWesixrcNnvjwBYlvX8H67pg0qYiFSAYzmFfLw8iY+X7SErv5i+kSHc3a8ZHRoHQmkp7PoNVrwNexaBqxdEXwGxt0L9tnZHF5FqxK4SdgUwxLKsMWW/Xw90tizr7j9c0weYjmOkLBlHIdt2qvdVCRORs5GZX8TUFXv5YEkix3KLiG0cwG29mjCgZR2cnQwc3uJY9HXLdCjOgwYdHGUsaiS4etodX0SqOLtK2JXA4L+UsE6WZd3zh2v8gFLLsrKNMcOANyzLan6C9xoLjAVo1KhRh71791ZIZhGpvnIKipm2Zj8fLdvDgWN5hAd7c0uPcK5oH4qnm7NjrtimrxwLv6bFg0ctaHsdxN4Cwc3sji8iVVSlvR15gtckAbGWZaWd7BqNhInI+SguKeXXbYd5f3Eimw4cJ8DLleu7hnFD18YE+7iDZUHSEljzIcT9DKXFENYT2l0PLS8BNy+7P4KIVCF2lTAXHBPz+wMHcUzMH/XH243GmLrAEcuyLGNMJ+BboLF1ilAqYSJSHizLYk3SMd5bnMi8uCO4OjtxSUx9bujamDYNazkuyjoCGz6F9VMhYy+4+0HryxyLwIZ21DIXInJadi5RMQx4HccSFR9ZlvWcMeYOAMuy3jXG3A2MA4qBPOB+y7KWn+o9VcJEpLztTs3m42V7+G79QXILS2jTsBY3dGnMRTH18HB1dkzk37sMNn4O23+EolwIjoC2oyDmGm2PJCInpcVaRUTOQFZ+Ed+tP8gnK5JITM0h0NuNqzs25LrOjQgNKLsNWZAF276HDZ/D/pWOjcOb9nesyh85DNx9bP0MIlK5qISJiJwFy7JYvjudT5YnMXfHEQD6tajDNR0b0icyBBdnJ8eFabsco2Obv4bMA47NwyOHQvSV0Kw/uLjb+ClEpDJQCRMROUcHM/L4fOVevl57gLTsAmr7unNlbChXxTakcZC346LSUti/CrZ84xglyzsKHv7QaoSjkDXu7thQXERqHJUwEZHzVFRSyvy4FL5es58FO1MotaBb0yCu7tiQwa3rOuaOAZQUQeJCRyGLmwmF2eBT11HIWg137FupQiZSY6iEiYiUo8PH8/l23X6mrd3P/qN5+Hu6cmnb+oxsH0pMqD/mf09NFuZC/K+wdTok/AYlBeAdAi0udhSysJ7au1KkmlMJExGpAKWlFisS0/lqzX5mbztMYXEpTUK8GdmuAZe2a/D/k/kBCrIhYQ7s+Ani50BRjmNB2BYXQcvh0LSv5pCJVEMqYSIiFex4XhG/bDnEd+sPsjrpKACdwwMZ2b4BQ6Pr4efxhxGvojzYPR+2/wQ7f4GC4+DmCxGDHRP7mw0Az1o2fRIRKU8qYSIiF9D+o7n8sOEg3284SGJaDu4uTgxoVYeR7RrQKyIE1/89XQlQXAh7FsOOHyFuFuSmgZOLY+5Y5FCIGAJBTe37MCJyXlTCRERsYFkWG/dn8P2Gg8zYlMyx3CJqebkypHVdLo6pT5cmgf+/3AVAaQkcXOcYHYv/FVK2O44HRzhGySKGQsPO4OxizwcSkbOmEiYiYrPC4lIWxafy8+Zk5m4/Qk5hCUHebgyJchSyTuGBODv9ZRukY0kQPxt2zoKkZVBaBJ4BjtuVTftD037gW8eWzyMiZ0YlTESkEskvKmHhzhRmbD7E/B0p5BWVEOLrzrCoulzcpj4dGgXg9NdClp/pmEcW/6vjScvcNMfxOtHQrJ+jlDXqosn9IpWMSpiISCWVW1jM/LgUZm4+xPy4FAqKS6nr58Gg1nUY3LouncID/zyHDByLwx7eDLvnwa75ju2TSovB1RvCejhW62/a3zGXTJuMi9hKJUxEpArILihm3o4jzNpyiEXxqeQXleLv6Ur/FrUZ1LouvSNC8HQ7wUKvBVmwZ0lZKZsHx/Y4jtdqBGG9ILwXhPcEv/oX9gOJiEqYiEhVk1dYwuKEVGZvO8y8HSkczyvCw9WJXs1DGNS6LgNa1qaWl9uJX3w00VHGEhdC0lLIz3AcD2zqKGPhvRwLxfrUvmCfR6SmUgkTEanCikpKWbPnKLO3HWbO9iMcOp6Ps5Ohc3gg/VvWoW9kCE1CfE784tJSOLLFMVKWtAT2LoeCTMe5kBaOMhbeExp1A5+QC/ehRGoIlTARkWrCsiy2HDzuKGTbjpCQkg1AWJAXfVvUpm9kbTo3CcTd5ST7U5YUw6FNkLTYUcz2rYCiXMe5wKaO9ckadXF815wykfOmEiYiUk3tP5rLgp0pLIhLYfnudAqKS/Fyc6Z7s2D6Rtamb4sQ6vl7nvwNigsheYNjcv++sq88x4r/eAWXFbKyUlY3BlxOcgtURE5IJUxEpAbIKyxhRWIa8+NSWBCXysGMPABa1vOjV0QwPZuFEBsWgIfrSUbJACwL0hIcI2T7Vjq+/2+iv4sH1GsLDTpAg/aO7wFhGi0TOQWVMBGRGsayLBJSsssKWQrr9h6juNTC3cWJTuGB9GgWTPdmwbSq5/f3Ncn+KuvI/4+UHVznuJ1ZnO845xVUVsrKvuq3B++giv+AIlWESpiISA2XU1DMqj3pLElIY2lC2u9zyYK83ejWLJiezYLp0TyY+rVOcevyf0qKHFsqHVwHB9Y5vqfGAWV/nwSEOcpYvRjHLcx6bcA7uMI+m0hlphImIiJ/ciQzn6UJaSzd5fhKzSoAIDzYm87hgXRuEkjn8KAzK2XgWKsseaOjkB1c6/j5+P7/P+9b/w+lrOx7rUa6lSnVnkqYiIiclGVZ7DySxdKENFYmprN6z1Ey84sBaBjoSefwIDqHB9KlSRANA73O/I1zjzpW9j+8BQ5tdvycFg9WqeO8h7+jjNWNgdotoXYrCIkE95MstyFSBamEiYjIGSsptdhxKJNVe46yKjGd1UlHycgtAqBBLU86hwfSKTyQ9o0DaBbic/o5ZX9UmOu4lXl48/8XsyP/1979x9Z13nUcf399f//0dfwjce00cZvQpVu6/phQoaJMFIkhJro/VjGgpZqQQKjAhkCwIhASf/EHIEAaY2gbdFq1MUonpknAoKCiiZGOpl3TNumSxmnsxMmNf9yf9rWv7Yc/nuPcexMnjVPnHsf5vKSj85zH514f+6tcf3Kec57zRusaM/BnyAYPBMEsWAZ+CGLXeFZOZAtRCBMRkeu2uur4QbHKoZOzHBqf4dDJWWbqSwDkklHuu72PB27v4/49Be7dXSCXjG3wG6zA3CkoHoULR/26eMyfNVv14Q/rgb4xH8gG74L+/TCwH/r3QaqwuT+wyCZSCBMRkU3jnGN8us7h0yVefmeOV07P8db5Ks75S7zu2pnzwWxPH/ffXmBsIINdz7VfK03/CKbimz6UFd/0AW32JLiV1n7pgVYg698XtPf7GwQ0r5mETCFMRERuqEqjyfcnShx+p8TLp30wqwbXleWTUe4ZLXBwtJd7Rno5ONrLSCF1fcEMfDibO+XnM5s5DjMnYPqEb9cvtPazCPTt8cGsbwx2jPlg1jfm+zW8KV2gECYiIl21uuo4caHG4XfmeO1MmSOTZY6dq9Bc8X9zdmTiHBzp5Z7R3mBdYGc+cf3BbM1CCWbebgtnx/323ClYqnbumxtuhbKOgLbXT6mhOzdlEyiEiYhI6BrNFd46Vw1CWYnXJsscL9ZYWfV/hwZzCQ4M5zmwK8eB4TzvG85x52CWWKTnvX9z5/zdmnPjMDvuQ9lcsJ4dh+rZzv3jWR/GendDYTf0jvr22nZmCHo24bhk21MIExGRLWlhaYU3pyocmSxx5EyFo1MVThRrLK34aSxiEWPfUI4DwzkO7MpfDGcD2cTmHkhzAUqnLw9opQkoT8JiuXP/SBzyI0FAW1tGW9v5EYglN/cY5aakECYiIjeN5soqJy/UOTpV4ei5CkenqhybqlAMJpQFGMgm2DeUYf9Qjv07s+wbzLJvZ5bB7CYMaa6nUfZhrDThJ6EtT7QCWnkCque4+MSANekByA/7iWrzw374MzcM+dta61Sfhj23OYUwERG56c3UFjl2rsrRqQo/OF/leLHGifM1qovLF/fpTcXYN5Rl/1CWfcGyf2eO4XxyY/OZbdTyElTOtEJZacIPcVamWuv56ctfF01CbldnUFsLadmdkB3ySyKvsHaTUggTEZFtyTlHsbrI8fM1ThR9MDterPF2sXZxLjOAZKyHvf0Z9vZn2DOQZqw/w96BDGMDGYZyN+js2aWWF/0Zs+oUVM4G7bWgttY31Tlx7ZpIojOUZYf8dWkXt3dCZtCv9cSBLUUhTEREbjmz9SVOFGscL1YZv1Dn1Eyd8ek6E7MLF685A0jHI+zpzzA2kPZBLQhne/rTN25480qcg4U5H8ZqRT/lRu28b9eKvr3WV5/msiFQgFi6FdIyg5De4e/2TPf7IdLMgO9LB33xjM6y3UAKYSIiIoGVVcfZ0gLj061gdmq6zqmZeSZm51lebf1dTER7GO1LsXtH2q/70h3tQjrW3ZDW8YMsw/wM1IuXBLViq68+4/eZn4bV5fXfJ5oMAlqwZAZaAS3T3wpv6X4f3pIFTYK7AVcLYdFuH4yIiEiYIj3G7h0+TD3MYMfXlldWmZxbYHymzsSsD2WTcwtMzM3zyukS5YVmx/7ZRJTRvhSjfemLYW2kkGRXb4rbepMMZBM37lq0SBRyb1EPBwAACmdJREFUO/3Cwavv65y/uWB+prXUp1sBbX62tT037rcXK1d+v3jW31SQKgTrtiW5Tt/aEkvprFsbhTAREZFANNLD3gE/JLmeSqPJ5KwPZZNzC0FI88t3356mvrTS+X49xs58ktsKSYZ7Uwz3Jv1SSHFbb4pdvUn6M/Ebe9MA+OCTKvil/85re83y4uWBbWHOT4i7MNe5FI9Bo+TD22rzyu8ZSawTzgqQ7PVLIh+085dsB+3I9ootGo4UERHZBM455uabnC0tcK7cYKq8wNlyg3PlBmdLC0wF7fbr0QDikR529SbZ1ZtkKJdgKJdkKJ+4rN2bCnHo81o5B835y0PauktbmGtULn+iwXpimSuEtHVC28XtvF8ncv4MXpcn2dVwpIiIyA1mZuzIxNmRifOBkd5193HOMVNfYqrU4GzZh7Wz5QWmSj6gvX6mTLFaZP6SM2oA8WhPEMwuD2qDQXswm6AvE9+cpwxcDzN/oX884yev3YjVFT8E2ij7UNYoX2G71Nqen4bZt1vbVzsLtyae84EsmYeDj8HDv3N9P+smUAgTERHpEjNjIJtgIJvg4Oj6QQ2gtrhMsdKgWF30y1o7WJ+4UON/3p6m0lj/YvtCOkZ/Jk5/NsFANk5/JsGOTNy3s4mOr+WTsRs/HHoteiKtIcrr4Zyf3qMjtJWDdQ0Wq8FSCZaqv9EgRAphIiIiW0w2ESU7mOWOwavP+dVornChukix2qBYWWS6vsRMbZGZ2hIz9UWma0u8da7KTH2G0vz6Z4miPf4MXiucxelLB0smRiEdpy8doy8dpxCs0/HI1hsaNfMX/sdSfgLcm4BCmIiIyE0qGYtcvNPz3TRXVpmbX/IBrS2kzdQWma0v+XZ9kXdO15mrN6ktXmFKC/zQaCHVCmY7MvF1w9paiCukYuRTsfCGSbcohTAREZFbQCzS468ly13bg8WXllcpLSxRmm8yV19ibr5Jab593WofL9Yozft92+dZu1Q6HqE3FSOfjPl1Kkq+YztYJ6Od26kYma149u09UggTERGRy/gbAa49tIG/8aC6uEyp3gxCmg9mpfklKo1lygtNKgtNv240OVNqcHSqSmWh2fEM0PVEeuzycJb0QS6XjJFLRMkmg3YySi7Ravv+KIlo5L3+WjaVQpiIiIhsCjPzwSgZ4/b+dx8ibbey6qg2mlQWli+GtPbA5gPccsf22dIC5YVlaotNGs3Vd/0e8WhPEM58QPvoPcP86o9f47xpN4BCmIiIiIQu0mP++rH09T0SqbmySq2xTLWxTKXhr2mrNpaptrUrjebFfaqNJolouNeoKYSJiIjITS8W6aEvE6cvc/M811K3KYiIiIiEQCFMREREJAQKYSIiIiIhUAgTERERCYFCmIiIiEgIFMJEREREQqAQJiIiIhIChTARERGRECiEiYiIiIRAIUxEREQkBAphIiIiIiFQCBMREREJgUKYiIiISAjMORf2MWyImV0A3unCtxoAprvwfeTaqSZbk+qy9agmW5PqsvV0oyZ7nHOD633hpgth3WJm/+ec+1DYxyEtqsnWpLpsParJ1qS6bD1h10TDkSIiIiIhUAgTERERCYFC2JX9bdgHIJdRTbYm1WXrUU22JtVl6wm1JromTERERCQEOhMmIiIiEgKFsEuY2UfM7C0zO2Fmnwn7eG4VZrbbzP7LzI6a2Rtm9qmgf4eZ/buZHQ/WfW2veTqo01tm9lPhHf32Z2YRM3vFzL4VbKsuITKzgpk9Z2bHgn8zP6KahM/Mfiv4/HrdzL5qZknVpfvM7EtmVjSz19v6NlwHM3vAzI4EX/srM7PNPlaFsDZmFgE+C/w0cDfw82Z2d7hHdctYBn7bOXcAeBB4KvjdfwZ4wTm3H3gh2Cb42ieA9wMfAf46qJ/cGJ8CjrZtqy7h+kvgX51z7wM+iK+NahIiMxsBfhP4kHPuA0AE/3tXXbrv7/G/03bXU4fPAb8C7A+WS9/zPVMI6/TDwAnn3Enn3BLwNeDRkI/pluCcm3LOHQ7aVfwflRH87/+ZYLdngI8F7UeBrznnFp1z48AJfP1kk5nZKPAzwBfaulWXkJhZHngY+CKAc27JOVdCNdkKokDKzKJAGjiL6tJ1zrn/BmYv6d5QHcxsGMg7577r/MXzX257zaZRCOs0Aky0bU8GfdJFZrYXuA84BOx0zk2BD2rAULCbatU9fwH8LrDa1qe6hOcO4ALwd8EQ8RfMLINqEirn3BngT4HTwBRQds59G9Vlq9hoHUaC9qX9m0ohrNN64726fbSLzCwL/BPwaedc5Wq7rtOnWm0yM/soUHTOvXytL1mnT3XZXFHgfuBzzrn7gDrB0MoVqCZdEFxj9CgwBtwGZMzs8au9ZJ0+1aX7rlSHrtRHIazTJLC7bXsUfzpZusDMYvgA9qxz7vmg+3xwWphgXQz6VavueAj4WTM7hR+e/wkz+wqqS5gmgUnn3KFg+zl8KFNNwvWTwLhz7oJzrgk8D/woqstWsdE6TAbtS/s3lUJYp+8B+81szMzi+Iv1vhnyMd0SgrtOvggcdc79eduXvgk8GbSfBP65rf8TZpYwszH8RZMvdet4bxXOuaedc6POub34fw//6Zx7HNUlNM65c8CEmd0VdD0CvIlqErbTwINmlg4+zx7BX9uqumwNG6pDMGRZNbMHg3r+UttrNk10s9/wZuacWzazXwf+DX9ny5ecc2+EfFi3ioeAJ4AjZvZq0Pf7wJ8AXzezX8Z/yD0G4Jx7w8y+jv/jsww85Zxb6f5h37JUl3D9BvBs8J/Fk8An8f+pVk1C4pw7ZGbPAYfxv+dX8LOxZ1FdusrMvgp8GBgws0ngj7i+z6xfw99pmQL+JVg291g1Y76IiIhI92k4UkRERCQECmEiIiIiIVAIExEREQmBQpiIiIhICBTCREREREKgECYico3M7MNm9q2wj0NEtgeFMBEREZEQKISJyLZjZo+b2Utm9qqZfd7MImZWM7M/M7PDZvaCmQ0G+95rZv9rZq+Z2TeCZwBiZvvM7D/M7PvBa+4M3j5rZs+Z2TEzezaYTVtEZMMUwkRkWzGzA8DPAQ855+4FVoBfBDLAYefc/cCL+Fm0Ab4M/J5z7h7gSFv/s8BnnXMfxD8DcCrovw/4NHA3cAf+aQ8iIhumxxaJyHbzCPAA8L3gJFUK/7DeVeAfgn2+AjxvZr1AwTn3YtD/DPCPZpYDRpxz3wBwzjUAgvd7yTk3GWy/CuwFvnPjfywR2W4UwkRkuzHgGefc0x2dZn94yX5Xe2bb1YYYF9vaK+hzVESuk4YjRWS7eQH4uJkNAZjZDjPbg/+8+3iwzy8A33HOlYE5M/uxoP8J4EXnXAWYNLOPBe+RMLN0V38KEdn29D84EdlWnHNvmtkfAN82sx6gCTwF1IH3m9nLQBl/3RjAk8DfBCHrJPDJoP8J4PNm9sfBezzWxR9DRG4B5tzVzsiLiGwPZlZzzmXDPg4RkTUajhQREREJgc6EiYiIiIRAZ8JEREREQqAQJiIiIhIChTARERGRECiEiYiIiIRAIUxEREQkBAphIiIiIiH4f7ranH4uXUFsAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 720x432 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plot_loss_curve(training_loss_list, testing_loss_list)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 更换数据集"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "我们换一个数据集，使用MNIST手写数字数据集。\n",
    "\n",
    "MNIST是最有名的手写数字数据集之一，主页：http://yann.lecun.com/exdb/mnist/\n",
    "\n",
    "MNIST手写数字数据集有60000个样本组成的训练集，10000个样本组成的测试集，是NIST的子集。数字的尺寸都是归一化后的，且都在图像的中央。可以从上方的主页下载。\n",
    "\n",
    "我们使用的数据集是kaggle手写数字识别比赛中的训练集。数据集一共42000行，785列，其中第1列是标记，第2列到第785列是图像从左上角到右下角的像素值。图像大小为28×28像素，单通道的灰度图像。\n",
    "\n",
    "我们使用的是kaggle提供的MNIST手写数字识别比赛的训练集。这个数据集还是手写数字的图片，只不过像素变成了 $28 \\times 28$，图片的尺寸变大了，而且数据集的样本量也大了。我们取30%为测试集，70%为训练集。训练集样本数有29400个，测试集12600个。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "X.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "from sklearn.datasets import load_digits #不会导入数据\n",
    "digits=load_digits()\n",
    "X = digits.data[:, 1:]\n",
    "Y = digits.data[:, 0]\n",
    "\n",
    "trainX, testX, trainY, testY = train_test_split(X, Y, test_size = 0.3, random_state = 32)\n",
    "\n",
    "trainY_mat = np.zeros((len(trainY), 10))\n",
    "trainY_mat[np.arange(0, len(trainY), 1), trainY] = 1\n",
    "\n",
    "testY_mat = np.zeros((len(testY), 10))\n",
    "testY_mat[np.arange(0, len(testY), 1), testY] = 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "len(X)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "trainX.shape, trainY.shape, trainY_mat.shape, testX.shape, testY.shape, testY_mat.shape"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "绘制训练集前10个图像"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "_, figs = plt.subplots(1, 10, figsize=(8, 4))\n",
    "for f, img, lbl in zip(figs, trainX[:10], trainY[:10]):\n",
    "    f.imshow(img.reshape((28, 28)), cmap = 'gray')\n",
    "    f.set_title(lbl)\n",
    "    f.axes.get_xaxis().set_visible(False)\n",
    "    f.axes.get_yaxis().set_visible(False)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## test：请你使用kaggle MNIST数据集，根据下表设定各个超参数，计算测试集上的精度，绘制损失值变化曲线，填写下表\n",
    "\n",
    "任务流程：\n",
    "1. 对数据集进行标准化处理\n",
    "2. 设定学习率和迭代轮数进行训练\n",
    "3. 计算测试集精度\n",
    "4. 绘制曲线"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###### 双击此处填写\n",
    "\n",
    "精度保留4位小数；训练时间单位为秒，保留两位小数。\n",
    "\n",
    "隐藏层单元数 | 学习率 | 迭代轮数 | 测试集精度 | 训练时间(秒)\n",
    "-|-|-|-|-\n",
    "100 | 0.1 | 50 | 0.6703755215577191 | 0.3141636848449707\n",
    "100 | 0.1 | 100 | 0.827538247566064 | 0.5589699745178223 \n",
    "100 | 0.1 | 150 | 0.9179415855354659 | 0.8267562389373779\n",
    "100 | 0.1 | 500 | 0.9624478442280946 | 2.6545965671539307\n",
    "100 | 0.01 | 500 | 0.6675938803894298 | 6.239083290100098"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "# YOUR CODE HERE\n",
    "stand = StandardScaler()\n",
    "trainX_normalized = stand.fit_transform(trainX)\n",
    "testX_normalized = stand.fit_transform(testX)\n",
    "\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "start_time = time()\n",
    "\n",
    "h = 100\n",
    "K = 10\n",
    "parameters = initialize(h, K)\n",
    "training_loss_list, testing_loss_list = train(trainX_normalized, trainY_mat, testX_normalized, testY_mat, parameters, 50, 0.1)\n",
    "prediction = predict(testX_normalized, parameters)\n",
    "print('testing accuracy:', accuracy_score(prediction, testY))\n",
    "plot_loss_curve(training_loss_list, testing_loss_list)\n",
    "\n",
    "end_time = time()\n",
    "print('training time: %s s'%(end_time - start_time))\n",
    "\n",
    "start_time = time()\n",
    "\n",
    "h = 100\n",
    "K = 10\n",
    "parameters = initialize(h, K)\n",
    "training_loss_list, testing_loss_list = train(trainX_normalized, trainY_mat, testX_normalized, testY_mat, parameters, 100, 0.1)\n",
    "prediction = predict(testX_normalized, parameters)\n",
    "print('testing accuracy:', accuracy_score(prediction, testY))\n",
    "plot_loss_curve(training_loss_list, testing_loss_list)\n",
    "\n",
    "end_time = time()\n",
    "print('training time: %s s'%(end_time - start_time))\n",
    "start_time = time()\n",
    "\n",
    "h = 100\n",
    "K = 10\n",
    "parameters = initialize(h, K)\n",
    "training_loss_list, testing_loss_list = train(trainX_normalized, trainY_mat, testX_normalized, testY_mat, parameters, 150, 0.1)\n",
    "prediction = predict(testX_normalized, parameters)\n",
    "print('testing accuracy:', accuracy_score(prediction, testY))\n",
    "plot_loss_curve(training_loss_list, testing_loss_list)\n",
    "\n",
    "end_time = time()\n",
    "print('training time: %s s'%(end_time - start_time))\n",
    "start_time = time()\n",
    "\n",
    "h = 100\n",
    "K = 10\n",
    "parameters = initialize(h, K)\n",
    "training_loss_list, testing_loss_list = train(trainX_normalized, trainY_mat, testX_normalized, testY_mat, parameters, 500, 0.1)\n",
    "prediction = predict(testX_normalized, parameters)\n",
    "print('testing accuracy:', accuracy_score(prediction, testY))\n",
    "plot_loss_curve(training_loss_list, testing_loss_list)\n",
    "\n",
    "end_time = time()\n",
    "print('training time: %s s'%(end_time - start_time))\n",
    "\n",
    "h = 100\n",
    "K = 10\n",
    "parameters = initialize(h, K)\n",
    "training_loss_list, testing_loss_list = train(trainX_normalized, trainY_mat, testX_normalized, testY_mat, parameters, 500, 0.01)\n",
    "prediction = predict(testX_normalized, parameters)\n",
    "print('testing accuracy:', accuracy_score(prediction, testY))\n",
    "plot_loss_curve(training_loss_list, testing_loss_list)\n",
    "\n",
    "end_time = time()\n",
    "print('training time: %s s'%(end_time - start_time))"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
